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ABSTRACT 



This thesis is a continuation of previous work to apply 
modern multiparameter estimation techniques to the problem 
of estimating attrition rates for a large number of small 
inventory cells in manpower planning models used by the U.S. 
Marine Corps. The main advances involve the promising 
introduction of empirical Bayes (non-constant shrinkage) 
techniques, recognition of the non symmetric nature of the 
errors with a response to this, and some insight into all 
aggregation plans that should help provide greater stability 
for the estimation methods. In addition, the roles of some 
middle level methodological choices are explored. 
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THESIS DISCLAIMER 
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The reader is cautioned that computer programs 
developed in this research may not have been exercised for 
all cases of interest. While every effort has been made, 
within the time available, to ensure that the programs are 
free of computational and logic errors, they cannot be 
considered validated. Any application of these programs 
without additional verification is at the risk of the user. 
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INTRODUCTION 



I . 



A. PURPOSE 

In support of project OPUS (Officer Planning and Utility 
System), Headquarters, U.S. Marine Corps, in 1985, requested 
assistance in exploring new methods of generating manpower 
loss rates in order to improve upon the currently used 
method of force rate analysis. In response, students and 
faculty in the Operations Research department at the Naval 
Postgraduate School began exploring the use of modern 
multiparameter estimation techniques for this problem. 
Special emphasis is placed on the small call problem, i.e., 
categories of officer skill, grade and length of service 
which have low inventory figures. Historically, rate 
estimators for small cells are unstable and a large number 
of these cells exist. 

This paper builds upon previous work on this problem. 
The main advances are to append new measures of 
effectiveness as requested by the Navy Personnel Research 
and Development Center (NPRDC) , introduce a class of 
empirical Bayes estimators, and explore the effects of some 
middle level choices in applying the new and existing 
methods. The ultimate goal is to refine the techniques 
presented here in order to validate a clear policy for 
predicting loss rates. 
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B . BACKGROUND 



For an introduction into Marine Corps policy concerning 
manpower planning, the reader is referred to a thesis 
submitted in September 1985, at the Naval Postgraduate 
School, by Major D.D. Tucker [Ref. 1], who presents a 
detailed background into the Marine Corps officer structure 
and the manpower planning process. Tucker also provides 
basic attrition rate theory and calculates overall attrition 
rates in several different formats as they pertain to the 
Marine Corps . As the aggregation rate begins to grow 
smaller with refinement. Tucker illustrates the irregular 
behavior of losses due to voluntary attrition and the small 
numbers of losses obtained in such aggregates. Tucker [Ref. 
l:p. 50] introduces the James-Stein technique of loss rate 
estimation as a shrinking of individual cell averages toward 
a grand average in order to reduce the risk, or squared 
differences of forecast and actual values, the goal being an 
improvement over the classical maximum likelihood estimator, 
referred to from now on as the MLE. Tucker chose the ranks 
of First Lieutenant and Lieutenant Colonel for evaluation 
using military occupational specialty (MOS) groups of combat 
support and ground combat (three MOS’s were selected for 
each group) . These aggregation schemes were carried through 
other studies and are included for historical purposes in 
this paper. 
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Another thesis subraicced at the Maval Postgraduate 
School in March 1986, by Major J.R. Robinson [Ref. 2] 
focused on a technique called limited translation of the 
James-Stein estimator. This particular model attempts to 
minimize the risk mentioned above by reducing the shrinkage 
of rate estimates for cells which fall outside a certain 
range of values centered on the grand mean. It was thought 
that shrinking all cells by the same rate toward this grand 
mean may be unwise since the attrition rates for those 
particular cells, those farthest from the grand mean, may 
contain attributes which are different than the majority of 
the cells in the aggregate. Additionally, Robinson’s work 
carries forward the study done by Tucker concerning MLE and 
James-Stein estimators. 

These feasibility studies were limited in scope due to 
the format of the summary date tape available (acquired from 
NPRDC), which allows only coarse cell definition. A cell is 
defined to be a cross-classification of the forty military 
occupational fields (OF), thirty-one lengths of service 
(LOS) , and ten grades for a total of 12,400 categories for 
manpower planning purposes [Ref. 2:p. 10]. The data does 
not distinguish between limited duty officers, those 
officers specifically designated for limited duty within 
certain MOS's, and unrestricted officers (regular and 
reserve) [Ref. l:p. 20]. These factors dictated a broad 
aggregation scheme that created stable cell inventories and 
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allowed the common variance assumption crucial to their 
work. But some of the application models require more 
refined information along those lines and experiments with 
finer aggregation levels, necessary in some real cases, 
caused this assumption to be violated. 

Tucker and Robinson did conclude that the current scheme 
of attrition rate prediction could be improved upon by their 
methods, although no dominant scheme was uncovered. Tucker 
and Robinson did identify problem areas in their studies 
that required additional work. Tucker [Ref. l:p. 71] cited 
a need for a better aggregation method to produce a more 
homogeneous attrition behavior in all cells. Also, small 
probabilities of loss within a cell were not dealt with 
successfully. Robinson [Ref. 2:p. 32] revealed the small 
cell problem and the inability to normalize the cell means 
or stabilize the variance with a data transformation. The 
recommendations of the previous studies highlighted the need 
for a more refined data tape with current information to 
include, among others, full Military Occupational Specialty 
(MOS) information, grade separation to include 
regular /reserve status, promotion zone data and breakdown of 
attritions by type. 

C. PROGRESS 

The project of loss rate estimation, sponsored by NPRDC , 
is currently moving into the implementation stage. An 
operational data tape is now available which includes 
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detailed information on Marine Corps officers for the years 
1977-86. Work has begun at the Naval Postgraduate School to 
break out the data in a useful format and its availability 
in now nearing completion. 

In addition to the results stated earlier, Tucker and 
Robinson provided loss rate estimates for small cells with 
no attrition, i.e., MLE equal to zero, and risk effects of 
the different schemes studied. These results provided the 
springboard for this study and the requirements for 
additional future studies. 

Captain R.W. Larsen [Ref. 3] has developed a promising 
aggregation scheme based on cluster analysis, i.e,, a 
classification scheme to aggregate cells which reflect a 
greater degree of homogeneity in attrition rates than are 

I 

allowed using current aggregation methods. Table 1 displays 
the resulting scheme, given in year of current service (YCS) 
groups proposed for implementation by Larsen. 



Table 1. AGGREGATION METHOD PROPOSED BY LARSEN 



MOS Catesory 


Bounded YCS Groups 


Fixed-Wing Pilots 


(1-6,8-19) 


(7) 


(20-25) 


(26) 


Rotary- Wing Pilots 


(1-5,8-19) 


(6,7) 


(20-25) 


(26) 


Naval Flight Officers 


(1-5,8-19) 


(6,7) 


(20-25) 


(26) 


La\v7ers 


(1-6,8-19) 


(7) 


(20-25) 


(26) 


All Else 


(1-3,6-19) 


(4,5) 


(20-25) 


(26) 
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The parentheses encompass YCS groups which behave similarly. 
In each case the grade is fixed, i.e., one should aggregate 
over YCS before aggregating over grade. 

This information will allow for more successful 
application of the loss rate estimation techniques of 
interest here when the new data tape is provided. 



D. GOALS 

In an attempt to refine the methods of attrition rate 
estimation presented in the previous pilot studiee, the 
following goals were set for this paper: 

1. Modify the existing methods to: 

a) use the refined data format, 

b) extend the work to the case of unequal or non 
constant variance in the small cell inventory 
problem , 

c) study alternate transformation inversion 
techniques , 

d) introduce and evaluate additional measures of 
effectiveness suggested by NPRDC to include, but 
not limited to, cross validation underages and 
overages of the deviation between forecasts and 
actuals , 

e) examine the graphical effect of shrinkage in the 
original scale. 

2. Introduce empirical Bayes estimation method in several 
forms for attrition rates for consideration as another 
option to solve the problem. 
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II. EXTENSIONS TO ESTIMATION METHODS 



A. GENERAL 

As stated in Chapter I, Tucker and Robinson published 
pilot studies introducing the James-Stein, limited 
translation James-Stein, Maximum Likelihood and Transform 
Scale Cell Average estimators as alternatives to the 
aggregate methods previously used for manpower planning in 
the Marine Corps. Their performance was promising, and this 
chapter will investigate ways to enhance their performance 
and possibly develop a dominant estimator for future 
attrition prediction. Additionally, NPRDC has expressed 
interest in additional measures of effectiveness for these 
estimators which may lead to a sharper direction for the 
goal of producing a viable estimation policy. This chapter 
concentrates on the following issues concerning previous 
studies . 



B. NON CONSTANT VARIANCE 

For a cell having an inventory n, and attrition rate p, 
the number of attritions y is described by the Binomial 
(n,p) distribution. The variance of the estimator of 
p = y/n is given by the familiar formula 



var(p) =Eli_:_Ei 
n 



. 1 ) 
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The Freeman-Tukey double arcsine transformation was used 
by Tucker [Ref. l:p. 55] and Robinson [Ref. 2;pp. 74-75] in 
order to map the raw losses to the transformed scale. This 
formula is given by 



= 0.5[n + 0.5] 1 ^2 (sin'i [ 



(n + i: 



- 1 ] 






( 2 . 2 ) 



Note that n=n(t) may change with time and thus, so would x. 
The question of averaging over time before or after applying 
the Freeman-Tukey transformation is not discussed in this 
paper. One of the goals of using this transformation is to 
stabilize the variance of this estimate. In the work of 
Tucker and Robinson, the assumption of normally distributed 
random variables with common variance was the setting in 
which the James-Stein estimator was expected to perform 
well. Indeed, for n moderate or large and p not too 
extreme, the variance of the Freeman-Tukey transformation is 
approximately one. But as discovered by Robinson, this 
assumption does not hold well when applied to cells with 
small inventory. A more careful look at the variance of the 
transformed data is taken here, with particular interest in 
the region of unstable variance where n and p are small. 

Some exploratory graphical work has revealed that the 
natural log of the variance of x, the transformed value. 
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behaves quite well as a linear function of u ( =E[X]) and la- 
1 in this unstable zone. This leads to the interpolatory 
formula : 



var(x) = 1.6835U-- - 1 ) • ^ b i , (2.3) 

1 . 001 ^ u ^ 2 . 2 

The upper limit, p = 2.2, is, the value for which var(x) = 1. 
Since the least squares fitting process was applied to the 
logarithms, p can not be allowed to fall to one or below. 
The lower limit, p = 1.001, is an arbitrary value that meets 
this constraint. The coefficients were validated by a 
linear regression on the natural log of this formula. The 
actual variance function for n = 1,2,3,4,5,7,10 appear in 

Figure 1, and the fit is remarkably good for n > 3. The 
upper graph shows the region of special interest for this 
paper and this portion of the curve is fitted by the formula 
in Equation 2.3 above. The lower plot exhibits the outer 
tale effect yet to be explored. The small values of n(<3) 
do not maintain the stable variance region, i.e., var[x]=l, 
for long before falling back. As the value for n increases, 
the more stable the variance becomes. Table 2 is a partial 
display of the residuals for selected n values, fifteen 
equally spaced values in the range 1.001 < p < 2.2, computed 
by the linear regression formula above. For values of n > 2 
the difference between actual and fitted values are 
considered acceptable for this work. 
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Figure 1. variance of Freeman-Tukey Transform as 
Function of It’s Mean for Selected n. 



Table 2. RESIDUAL LISTING FOR SELECTED VALUES OF N 



N= 1 


N = 2 


N = 4 


N = 5 


N= 10 


% 


0.07191 


0.03912 


0.02139 


0.017 


0.008313 




0.08466 


0.02569 


0.002853 


-0.001695 


-0.01106 




0.1033 


0.0239 


-0.005448 


-0.01049 


-0.0206 




0.122 


0.03002 


-0.00655 


-0.01185 


-0.02257 




0.1385 


0.04068 


-0.004026 


-0.009566 


-0.02067 




0.15 


0.05388 


0.0000676 


-0.005727 . 


-0.0166 




0.1552 


0.06853 


0.004027 


-0.002103 


-0.01291 




0.1527 


0.08234 


0.007447 


0.0009541 


-0.009981 




0.1414 


0.09481 


0.00948 


0.002557 


-0.008299 




0.1209 


0.1055 


0.009825 


0.002448 


-0.008128 




0.08958 


0.1129 


0.008369 


0.0004583 


-0.009929 




0.04832 


0.1169 


0.004997 


-0.003511 


-0.01327 




-0.005254 


0.1171 


0.0001116 


-0.00893 


-0.0182 




-0.06914 


0.113 


-0.006798 


-0.01654 


-0.02463 




-0.1466 


0.104 


-0.01534 


-0.02585 


-0.03366 





The boundaries for y, in the unstable region, as a 
function of n and p are provided in Figure 2. The formula 
derived was quite adequate for this study. 

C. ALTERNATE TRANSFORM INVERSION METHODS 

There are three issues in the inversion of the data from 
the transformed space. The first is that of selecting 
values of n, i.e., an average of the n(t)’s appearing in the 
estimation year inventories. This issue is discussed later 
in the chapter, with the development of a simulation for 
comparing mean inventory values. The second is that of the 
inversion of the sum of the two arcsine functions. The 
third involves the case of time averaging in the transformed 
scale and the question of what matching form of average 
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Figure 2. Envelope for the Non-Constant Portion of the 

Fitted Variance as Functions of Attrition Rate 



inventory should be used in the inversion. The third issue 
is not covered in this paper, and is still pending. 
Treatment of the second issue will be provided first. 

The method used in the pilot studies of Tucker and 
Robinson, and the most transparent transform inversion 
formula is: 

Pi = 0.5[1 - sin(^^^^ 1/2 ) 1 / 2 ^ ^ (2.4) 

This will be referred to as the Basic Inversion. The 
formula offered by Carter and Rolph [Ref. 4] is given by: 



]pi 



sin2 (01 ) + 



(1 - 2sin20 t ) 
(4(n/k) + 2) 



1 3i 



( 2 . 



5) 



where ©i is x/ (n + 1 / 2 ) and Bi is the estimated amount of 
shrinkage in the transformed scale depending on the Bayes 
model used (constant or proportional prior) , k is the number 
of cells and n is inventory. 

Apparently, Carter and Rolph invert prior to shrinking 
and the above formula provides a shrinkage in the original 
scale. Rewriting the above, it appears to shrink towards 
^ = 1/2, and is given by: 



= l/2( 



2Bi 



4(n/k) + 2 



) + sin2 (01 ) ( 1 - 



23i 



4(n/k) + 2 



) 



( 2 . 6 ) 



An exact Freeman-Tukey inversion formula developed by 
Miller [Ref. 5] was recently uncovered. It is exact in the 
sense that for a single cell, the formula given by: 
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p(tj = 0.5fl - sgn(cos t) [1 - (sin t + 

s_in — 5 — 2 — i lZ - Rj- . r* : — %.'! 2 ] I / 2 ] (2.7) 
n 

where t = x/(n + l/ 2 )i ^2 + n/2 and a is the average central 
inventory, returns the empirical rate y/n. This will be 
referred to as the Freeman-Tukey Exact (FTE) inversion. 
Since this formula is applied using an average value of n, 
there may be concern that it oscillates between successive 
integer values of n. Figure 3 indicates, using values of 
n=5 and 10, that it is quite smooth as an interpolation 
formula. In the article, Miller suggested using the 
harmonic mean of the n values provided by the time changes. 
This issue is examined later in this chapter. The Miller 
formula was incorporated in the Marine Corps Inventory 
Attrition Analysis Program provided in Appendix B. 

D. SHRINKAGE VISUALIZATION 

When shrinkage is a fixed amount for all cells, it is 
readily visualized, for all degrees of shrinkage, by the 
linear diagrams presented in the works of Tucker [Ref. l:p. 
52] and Casella [Ref. 6]. Robinson [Ref. 2:p. 20] also 
presented diagrams of the limited translation option and how 
the shrinkage is affected. All of these diagrams refer to 
the transformed scale. 

Of more immediate interest is the question of how to 
visualize shrinkage diagrams of this type in the original 
scale. Illustrations of diagrams of this type are provided 



14 





15 



Figure 3. Suitability of the Exact Freeman-Tukey Inversion 
as an Interpolatory Formula. 



using both the basic inversion formula used by both Tucker 
and Robinson, and the FTE formula. Figure 4 shows how the 



diagram by Casella 


[Ref. 


6 : p . 


84] 


looks 


in the original 


scale using the 


same 


value 


of 


n for 


all cells . An 



interesting note to consider is the effect of the difference 
between inversion formulas, i.e., in the graph which 
combines the Basic and FTE inverses, the FTE inversion tends 
to shrink at a slower rate, as values get further from the 
mean. For a modest transition into the case of non constant 
n. Figure 5 provides for two values of n. Finally, Figure 6 
illustrates the point with many different values of n, taken 
from, actual Marine Corps officer inventory data. In this 
case, note that Basic inversion shrinks to one original 
scale value while the FTE inverse is not as focused. 

E. CHOICE OF AVERAGE INVENTORY VALUES OVER TIME 

A Fortran simulation was developed, included as Appendix 
C, to test the performance of each of several common means 
(arithmetic, geometric and harmonic) in the comparison of 
the inverted estimates with the empirical y/n. Values of n 
were produced from a Poisson distribution with rate X (1 < X 
< 20), and attrition values, y, were produced from the 
Binomial (n,p) distribution where a was in the range .01 < p 
>.40. The n and y values were then transformed by using the 
Freeman-Tukey formula. An exhaustive study of inversion 
variations with both the FTE formula and the basic formula 
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BASIC INVERSE (N-10) 
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Figure 4. Non Linear Shrinkage Diagrams for a Single 
Value of Inventory, 



BASIC INVERSE (N= 10,20) 
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Figure 5. Non Linear Shrinkage Diagrams for Two Values 
of Inventory. 



BASIC INVERSE 
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Figure 6. Non Linear Shrinkage Diagrams with Empirical 
Bayes Estimators for Test Case. 



used by Robinson were compared. Table 3 contains selected 
outputs which are typical of the broader study and shows how 
the means performed against each other in three different 
computations listed below. 

1. The BASIC Method utilizes the Basic inversion formula 
and the average transformed scale value. 

2. The FTEl Method utilizes Equation 2.7 with the t 

values equal to the average (over replications, i.e., 
time) of x/(n + 1 / 2 ) given by Equation 2.2. 

More explicitly, see Equation A. 16. 

3. The FTE2 Method utilize Equation 2.7 with the t 
values equal to the average (over replication, i.e., 
time) of X directly, and x given by Equation 2.2 
followed with a division by (n + 1/2) 

In the context of this study, the arithmetic mean tended 

to perform as well or better than the geometric or harmonic 

means. The choice of arithmetic mean also circumvents the 

problem of an undefined term when ni = 0 is encountered 

while using the harmonic of geometric means. 

F. QUESTIONS RELATING TO VALIDATION 

The pilot studies by Tucker and Robinson both used four 
years (1977-80) for estimation of the attrition rates and 
the remaining three years (1981-83) for validation. The 
estimation rates were used in an attempt to forecast the 
values for the three validation years. The quality of the 
forecasts deteriorates as lead time increases in the future. 
This was anticipated, in fact, Rowe et al., at NPRDC have 
shown a time series effect in their work forecasting 
attritions in the Navy [Ref. 7] and are currently 
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Table 3. SIMULATION OUTPUT OF MEAN VALUES FOR INVERSION 
UTILIZATION 





REPS 


P 


METHOD 


ARITHMETIC 


GEOMETRIC 


HAR.MONIC 


2.0 


40 


0.05 


BASIC 


0.1068 


0.1176 


0.1308 








FTEI 


0.0252 


0.0189 


0.0122 








FTE2 


0.0181 


0.0208 


0.0245 


2.0 


40 


O.IO 


BASIC 


0.1399 


0.1537 


0.1708 








FTEI 


0.0535 


0.0458, 


0.0369 








FTE2 


‘ 0.0489 


0.0558 


0.0648 


2.0 


40 


0.20 


BASIC 


0.2402 


0.2629 


0.2906 








FTEI 


0.1900 


0.1829 


0.1739 








FTE2 


0.1657 


0.1869 


0.2144 


2.0 


40 


0.40 


BASIC 


0.3713 


0.4038 


0.4429 








FTEI 


0.3442 


0.3403 


0.3353 








FTE2 


0.3325 


0.3715 


0.4213 


4.0 


40 


0.05 


BASIC 


0.0750 


0.0815 


0.0921 . 








FTEI 


0.0235 


0.0200 


0.0147 








FTE2 


0.0200 


0.0221 


0.0258 


4.0 


40 


0.10 


BASIC 


0.1280 


0.1387 


0.1564 








FTEI 


0.0768 


0.0724 


0.0655 








FTE2 


0.0726 


0.0798 


0.0922 


4.0 


40 


0.20 


BASIC 


0.1965 


0.2125 


0.2386 








FTEI 


0.1446 


0.1406 


0.1339 








FTE2 


0.1485 


0.1628 


0.1869 


4.0 


40 


0.40 


BASIC 


0.4382 


0.4694 


0.5185 








FTEI 


0.4564 


0.4558 


0.4549 








FTE2 


0.4280 


0.4638 


0.5223 


8.0 


40 


0.05 


B.ASIC 


0.0552 


0.0580 


0.0616 








FTEI 


0.0181 


0.0166 


0.0148 








FTE2 


0.0156 


0.0165 


0.0177 


8.0 


40 


0.10 


BASIC 


0.1226 


0.1285 


0.1363 








FTEI 


0.0922 


0.0904 


0.0880 








FTE2 


0.0828 


0.0873 


0.0933 


8.0 


40 


0.20 


BASIC 


0.2035 


0.2131 


0.2255 








FTEI 


0.1778 


0.1762 


0.1741 








FTE2 


0.1708 


0.1797 


0.1915 


8.0 


40 


0.40 


BASIC 


0.3799 


0.3962 


0.4170 








FTEI 


0.3700 


0.3693 


0.3684 








FTE2 


0.3662 


0.3837 


0.4065 
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recommending a second order auto regression scheme for this 
purpose. In this regard, the work of Tucker and Robinson 
with validation may be viewed as forecasting hy persistence . 
and in light of the NPRDC work, can hardly be expected to 
perform well in an absolute sense. Their work does, 
however, provide some comparisons that can be useful in 
choosing among alternative techniques, and much can be 
learned from them. 

The first measure considered is the Mean Squared Error 
(MSE) of forecasts. The shrinkage estimation theory 
provides that these values should be about unity if 
persistence, independence over time, was agreeable. That 
is, useful validations can be obtained without modeling the 
time dependence of the attrition process. Only for a very 
few cases does this hold (i.e.. Warrant Officers, see 
Tucker). However, if one method of estimation consistently 
produces noticeably smaller MSE values than do the other 
methods, this would provide sufficient support for that 
method's continued development. 

The second measure is designed to consider performance 
in the original scale and is patterned after a chi-square 
statistic 



actuals - forecasts)^ 
forecasts 



If the k cells are independent, then 
approximately a chi-square random variable 



( 2 . 8 ) 

equation 2.8 is 
with k degrees of 
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freedom. Hence, the man value is about k and the standard 
deviation is (2k)^^^. Once again, we generally do not see 

I 

values this small, but use it for comparison purposes. As a 
side note, the value of k when using the MLE is generally 
smaller than that of the other estimators due to the 
omission of those cells with MLE equal to zero over all 
estimation years, i.e., zero forecasts. This provides and 
unnatural advantage for the MLE, which must be kept in mind 
when interpreting the results. 

Two disturbing features have emerged from the 
exploratory analysis. Often a technique that performs well 
for some cases using the MSE measure is not comparably 
supported by the chi-square measure. It seems the two 
measures do not weight the common features equally, and 
further study is necessary to develop understanding of the 
results. 

Discussion with personnel at NPRDC concerning these 
matters has led to the introduction of a third measure, 
which will be called the Mean Absolute Deviation (MAD). 
This value is computed in the original scale. Moreover, 
since the cost structure of over-forecasting is not 
symmetric, we separate the two parts into the average 
overage deviation given by: 

1 ^ 

ODEV = ^ SI [forecasts - actuals]* (2.9) 

^ 1 
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and average underage deviation 



UDEV = ^ 21 [actuals - forecasts]* (2.10) 

^ 1 

where the "+" refers to the positive part, as computed 
separately. Finally, 

MAD = ODEV + UDEV (2.11) 

This measure also has shown that in most cases, the 
overage exceeds the underage value, which in turn may lead 
to forecasting the need for an excess number of Marine 
officers . 
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III. EMPIRICAL BAYES ESTIMATION 



A . BACKGROUND 

Although the studies of Tucker and Robinson showed 
improvement over the current methods, no clear procedure was 
established and no plans for implementation were drafted. 
Problem areas were identified while working with small cell 
inventories and low loss rates. In the transformed scale, 
the James-Stein assumption of normal distribution of the 
observations with equal variance did not hold well in such 
cases. An alternate method of loss rate estimation is 
required to deal with this special case. 

B . INTRODUCTION 

Empirical Bayes is simply a generalization of the James- 
Stein estimator for use with small cell inventories and 
unstable variance. Random variables are assumed to be of 
the form: 



Xi ■ N(0i ,di ) , i=l,2 ,k (3.1) 

where ©i are the unknown parameters to be estimated and di 
are known variances. The Bayesian assumption is that the 
01 's are themselves from some distribution, in this case 
also from the normal distribution with mean, v, and 
variance, | 0 i a. The proportionality constant pi is 
normalized so that 
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s pi 



k, i=l,2,. 



• • t 



k 



(3.2) 



where k is the number of observations or cells. Thus, 

©1 ' N(v,^ia), i=l,2,...,k. (3.3) 

This distribution is known as the prior distribution of ©i 
and represents the additional assumed information on the 
value of © before observing X whose distribution depends on 
© [Ref. 8:p. 553]. The use of this information produces an 
estimate which is a weighted average of the p.ior mean, v, 
and the sample esti.mate x. In the empirical Bayes context, 
the values of v and a, are also estimated from the data 
since all the information necessary is contained in the 
marginal distribution of Xi (unconditional on ©t ) [Ref. 6:p. 
83] . This chapter focuses on that method. 

For a more detailed explanation of the Bayesian family 
of estimators, the reader is referred to Casella [Ref. 6], 
Carter and Rolph [Ref. 4] and Efron and Morris [Ref. 9]. 

C. THE ESTIMATION METHOD 

Fay and Harriot [Ref. 10] discuss the empirical Bayes 
technique and provide an estimate for the prior mean as 

H XI / (a* + dt ) 

V* = _] , i = l,2 k (3.4) 

x: 1/ (a* + di ) 

1 

where a* is an unbiased estimate of a. The empirical Bayes 
estimate of ©, commonly called the mean of the posterior 
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where a* is an unbiased estimate of a. The empirical Bayes 
estimate of ©, commonly called the mean of the posterior 
distribution becomes 



5E B = 



+ di 



) XI + ( 



di 



{a* + di ) 



) V* - i=l , 2 , . . . , k (3.5) 



This is the mean of the conditional distribution of © given 
X [Ref. 8:p. 556]. The weights of each of the prior and 
sample means are determined by the prior and observed 
variances . 

As a measure of the worth of this estimator, when the 
number of x's is > 4 [Ref. 11] , the risk, or sum of squared 
differences between the true unknown parameter and the 
empirical Bayes estimate, is less than that using the 
observed x's, for all ©i . Efron and Morris [Ref. 9] provide 
a rigorous proof. 

An important result of the empirical Bayes estimator is 
that shrinkage values are no longer constant, but depend 
upon information contained in each cell. Additionally, the 
point toward which the cell values are shrunk is the 
weighted mean, not the grand mean of the observations. 
Appendix A provides the derivation of the empirical Bayes 
estimator, variable shrinkage factors and calculation of 
the prior variance used in this study as they pertain to the 
Marine Corps officer attrition rates. 
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IV. VALIDATION PROCEDURES AND RESULTS 



A. GENERAL 

The outputs displayed in this chapter are excerpts from 
that produced by the Marine Corps Personnel Inventory 
Attrition Analysis (MCIAAP) (see Appendix B). Initially, 
comparisons will be made with the output displayed by 
Robinson to identify the validity of the Fortran code with 
the .iPL code used by Robinson, then new approaches will be 
displayed in an attempt to improve the results. The Fortran 
code in the MCIAAP incorporates most of the procedures 
studied by Tucker and Robinson. Specifically, the MLE , 
Transformed Scale Cell Average (TSCA) and the James-Stein 
(JS) estimators, shrinkage procedures and a basic inversion 
formula. The FTE inversion formula, the new measures of 
effectiveness from NF'RDC , are also included. For the MAD 
measure of effectiveness, only maximum likelihood and the 
modified empirical Bayes variables will be presented, since 
these showed the most promise over all. All other variables 
are displayed by the MCIAAP if interest is warranted by the 
reader. Chapter II describes all of these methods in some 
detail . 

B. HISTORICAL COMPARISONS 

Tables 4-9 were developed using the same aggregation 
scheme and format as Robinson [Ref. 2:pp. 39-41], although 
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not all of the variables are included (MLE, TSCA, and James 



Stein are the only variables duplicated in this paper) . 
Additionally, Robinson presented the number of cells, as k 
values, used for evaluation of each estimate in the 
aggregate Figure of Merit (FOM) display. These k values are 
not repeated in the current FOM tables, but are included, as 
parenthetical values, and discussed in the new MAD measure 
of effectiveness tables. Since this is the case here, the 
reader is reminded that the apparent competitive performance 
of the MLE in the following cases is due, in part, to the 
possible advantage of having a smaller number of cells used 
in this MOE ' s computation. 

Since the computations were identical for both studies 
for these tables, the output was expected to be nearly 
identical or some pattern discernable between the two. In 
very few cases was this found to be true, and after 
verifying that the algorithms coded both in the APL and 
Fortran versions were correct, we have no reasonable 
explanation for the discrepancies of the programs. It is 
suspected that the transfer of raw data into an APL 
workspace, accomplished by Tucker and Robinson, may be at 
variance with that into Fortran files, done by contractor. 

1. Aviation Aggregate Case 

In the case of aviation FOM ' s , displayed in Table 4, 
1981 values for LtCol's in the transformed scale are 
identical to Robinson’s [Ref. 2;p. 39] while the last two 
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Table 4. AVIATION FIGURES OF MERIT 





1981 1982 1983 

TRANSFORMED FOM 


« 


1st Lt 


MLE 


4.0863 


10.6590 


11.1338 




TSCA 


3.7420 


10.3636 


10.8711 




JS 


3.9694 


10.5671 


11.1626 




MOD JS 


4.7352 


11.2615 


12.1206 




EB 


4.0753 


10.7015 


11.3303 




MOD EB 


3.1833 


5.0568 


6.2870 




Li Col 


MLE 


4.3664 


8.5033 


8.1506 




TSCA 


5.7768 


10.2994 


10.3353 




JS 


5.7369 


10.2468 


10.2972 




MOD JS 


5.7208 


10.1865 


10.2756 




EB 


5.8482 


10.2298 


10.3534 




MOD EB 


4.3301 


3.9824 

ORIGINAL FOM 


2.7186 




1st Lt 


MLE 


21.0501 


43.1546 


53.0538 




TSCA 


56.4603 


93.7810 


98.4715 




TSCAM 


53.7884 


90.5164 


95.3329 




JS 


58.0150 


95.2150 


100.1664 




JSM 


55.3844 


91.9406 


97.0163 




MOD JS 


63.0345 


99.8419 


105.4897 




MOD JSM 


60.6413 


96.5615 


102.3172 




EB 


58.8824 


96.2449 


101. 2482 




EBM 


56.2111 


92.9699 


98.0958 




MOD EB 


19.7683 


38.4637 


49.7393 




MOD EBM 
Lt Col 


20.7987 


36.3673 


48.1659 




MLE 


33.9315 


59.4200 


22.1156 




TSCA 


38.4478 


57.0987 


69.0641 




TSCAM 


35.0390 


48.9374 


59.6957 




JS 


37.8248 


55.1913 


67.1696 




JSM 


33.2016 


46.9901 


56.7135 




MOD JS 


38.0291 


51.5325 


63.3908 




MOD JSM 


34.5466 


45.1554 


52.8285 




EB 


40.3421 


54.9389 


67.7091 




EBM 


36.7498 


46.9657 


57.8567 




MOD EB 


40.2056 


70.5381 


60.5749 




MOD EBM 


185.7476 


65.1553 


168.3976 
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years are smaller. The IstLt FOM ' s are larger for all three 
years. The modified empirical Bayes (Mod EB) values were 
best overall yet still may be too large (ideal values tend 
to be close to unity) . It is interesting to note that in 
most cases, the risk values tend to increase with future 
forecasts, yet the Mod EB FOM decreases as future forecasts 
increase with time. 

In the original scale, MLE tended to be smaller than 
Robinson's for all years, while TSCA and James-Stein (JS) 
were all larger. The FTE inversion worked best for all 
variables except for Mod EB in the first and third 
validation years. MLE and Mod EB worked best for IstLt 's 
and are approximately within the range of one standard 
deviation of the number of cells in the computation. For 
LtCol's, FTE inverted JS was best with MLE a competitive 
second, yet all were not in the desired range. 

In the MAD measure. Table 5, results are mixed 
between future forecasts (no set pattern) and variables for 
the two ranks (FTE Mod EB was best for IstLt 's and MLE was 
best LtCol's) . The mix in overage and underage is most 
desirable in the first year for all cases, yet the number of 
cells (k) forecasted as overage exceed the underage forecast 
cell numbers in most cases. This is cause for concern and 
the performance of these estimates in this aggregate is 
suspect. The second year values for LtCol estimates 
displays a preferred ratio of k values, but the total MOE 
value may still be too large. 
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Table 5. MEAN ABSOLUTE DEVIATION (AVIATION) 





1981 


1982 


1983 


1st Lt 
MLE 


UNDER.AGE 


0.1256( 3) 


0.0439( 1) 


0.3069( 1) 


OVERAGE 


2.3985( 8) 


3.9267(10) 


4.2522(10) 


TOTAL 


2.5240(11) 


3.9706(11) 


4.5591(11) 


MOD EB 
UNDERAGE 


0.1822( 4) 


0.04 10( 1) 


0.3281( 1) 


OVERAGE 


2.3625( 8) 


3.9040(11) 


4.2351(11) 


TOTAL 


2.5447(12) 


3.9451(12) 


4.5633(12) 


MOD EBM 


UNDERAGE 


0.25 1S( 5) 


0.0580( 2) 


0.3420( 1) 


OVER.AGE 


2.1654( 7) 


3.6427(10) 


3.9677(11) 


TOTAL 


2.4172(12) 


3.7007(12) 


4.3098(12) 


Lt Col 
MLE 


UNDERAGE 


0.6005( 5) 


1.5139(10) 


0.0396( 1) 


OVERAGE 


0.7925( 8) 


0.3928( 3) 


1.3467(12) 


TOTAL 


1.3930(13) 


1.9067(13) 


1.3863(13) 


MOD EB 


UNDERAGE 


0.6053( 6) 


1.4717( 9) 


0.0841( 2) 


OVER.AGE 


1.3083(10) 


1.0970( S) 


2.2331(14) 


TOTAL 


1.9136(16) 


2.5687(17) 


2.3173(16) 


MOD EBM 


UNDER.AGE 


0.6890( 6) 


1.5366(10) 


0.1167( 2) 


OVERAGE 


1.0120(10) 


0.7179( 7) 


1.7501(14) 


TOTAL 


1.7010(10) 


2.2546(17) 


1.8669(16) 



2. Combat Support Aggregate Case 

Table 6 displays the combat support aggregate output 
from the program. In the transformed scale comparison with 
Robinson [Ref. 2tp. 401 , the values for the two ranks take 

opposite directions. Where the IstLt FOM values are all 

. » 
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Table 6. COMBAT SUPPORT FIGURES OF MERIT 



1981 1982 , 1983 

TR.ANSFORMED FOM 



1st Lt 



MLE 


2.1877 


2.9028 


2.2447 


TSCA 


1.6329 


2.4580 


1.8696 


JS 


1.6925 


2.4816 


1.8987 


MOD JS 


2.3684 


3.0520 


2.5356 


EB 


2.3323 


3.1991 


2.6868 


MOD EB 


2.2293 


2.0048 


' 2.1531 


Lt Col 


MLE 


1.3103 


1.6978 


2.3108 


TSCA 


0.7502 


1.4070 


1.7540 


JS 


0.6421 


1.2345 


1.5429 


MOD JS 


0.7009 


1.0929 


1.3082 


EB 


0.7736 


1.2589 


1.5799 


MOD EB 


1.8691 


2.1236 


4.9017 



ORIGINAL FOM 

1st Lt 



MLE 


111.4260 


112.1917 


70.2767 


TSCA 


73.5008 


90.0523 


73.9380 


TSCA.M 


79.0947 


93.8244 


80.9764 


JS 


76.9945 


89.8064 


71.5890 


JS.M 


429.8503 


84.4175 


65.2639 


MOD JS 


143.1700 


113.0083 


90.6783 


.MOD JS.M 


168.9385 


114.7083 


92.6599 


EB 


102.9525 


1 19.0900 


96.3373 


EB.M 


108.3578 


111.2701 


86.8138 


MOD EB 


68.7251 


63.9054 


48.2719 


.MOD EB.M 


176.3961 


59.5760 


43.0037 


Lt Col 


MLE 


36.0152 


66.0261 


42.4097 


TSCA 


61.7769 


52.0346 


99.7917 


TSCAM 


65.9535 


49.6307 


107.1016 


JS 


61.0511 


49.2280 


96.7861 


JSM 


66.3509 


47.8490 


118.2564 


.MOD JS 


66.5470 


49.9366 


99.2566 


MOD JSM 


81.4722 


50.2637 


103.8460 


EB 


69.7508 


55.9812 


104.5311 


EBM 


71.5728 


50.0281 


109.2012 


MOD EB 


62.6316 


50.3125 


96.8693 


MOD EBM 


83.1337 


48.6576 


127.1541 



V t 
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larger (except for MLE in 1921) , the LtCol values are 
smaller. Overall, all estimates are competitive with each 
other over all validation years with TSCA performing best 
for the IstLt's and James-Stein outperforming the rest in 
the LtCol 's. Note here also that the Mod EB estimate did 
not fair well in either case, and was the worst estimate in 
risk for the LtCol 's. 

In the original scale , a different and more random 
effect is noticed between Table 6 and the results of 
R-obinson. Risk values tended to fluctuate with forecast 
lead time between estimates for IstLt's where the LtCol 
values were more stable, yet were quite different from the 
transformed scale pattern. In this case. Mod EB performed 
well for IstLt’s, and was the only one within the desired 
range of standard deviations. MLE for the LtCol ' s was also 
the only estimate desirable for this measure. The Miller 
inverse showed poorly in this aggregate for all estimation 
methods, in every case exceeding it's counterpart. 

In the MAD measure of effectiveness, Table 7, MLE 
performed the best in both ranks. The Mod EB effect is 
competitive in '•he IstLt's case, and all variables in the 
LtCol case are closely grouped. The number of cells are 
heavily weighted towards an overage forecast in all three 
validation years, throughout this aggregate. Once again, 
this is not a desirable quality for a MOE. 
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Table 7. MEAN ABSOLUTE DEVIATION (COMBAT SUPPORT) 





1981 


1982 


1983 


1st Lt 
MLE 


UNDER.AGE 


0.3737(14) 


0.1617(8) 


, 0.1781(8) • 


OVERAGE 


0.6731(19) 


0.9574(25) 


0.8835(25) 


TOTAL 


1.0468(33) 


1.1191(33) 


1.0616(33) 


MOD EB 


UNDER-AGE 


0.4028(15) 


0.1275( 9) 


0.1427( 9) 


OVERAGE 


0.7311(40) 


0.9972(45) 


0.9422(38) 


TOTAL 


1.1339(55) 


1.1246(54) 


1.0849(47) 


MOD EBM 


UXDER.AGE 


0.4592(16) 


0.1603( 9) 


0.1856(10) 


OVEFU\GE 


0.5930(39) 


0.8318(45) 


0.7575(37) 


TOTAL 


1.0522(55) 


0.9921(54) 


0.9431(47) 


Lt Col 
MLE 


UNDER.AGE 


0.1360( 9) 


0.1688( 7) 


0.054 1( 5) 


OVEFUXGE 


0.1510(14) 


0.2591(16) 


0.3415(18) 


TOTAL 


0.2870(23) 


0.4280(23) 


0.3957(23) 


MOD EB 


UNDERAGE 


0.1638(11) 


0.1721( 7) 


0.1162( 6) 


OVERAGE 


0.2232(29) 


0.3606(35) 


0.4548(29) 


TOTAL 


0.3870(40) 


0.5327(42) 


0.5710(35) 


MOD EBM 


UNDER-AGE 


0.2044(11) 


0.1901( 7) 


0.140S( 6) 


OVERAGE 


0.1263(29) 


0.2240(35) 


0.2589(29) 


TOTAL 


0.3307(40) 


0.4142(42) 


0.3996(35) 
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Ground Combac Aggregate Case 



In the final aggregate studied by Robinson [Ref. 
2:p. 41], that of ground combat MOS’s, the results again, as 
displayed in Table 8, were mixed between computational 
methods. In the original scale, value differences tended tc 
be larger than Robinson's in all cases for IstLt's, showing 
a proportional growth over time. Conversely, LtCol values 
were smaller in the first two validation years and slightly 
larger in the final year. TSC.-. looks good for the IstLt's 
with the others being competitive, yet the size of the 
values may not be desirable. Mod EB turned out the best for 
LtCol 's, but all other estimates are competitive and in the 
desired range. 

On the original scale, all values for the estimates 
for IstLt's are too large, with Mod EB performing well in 
that regard. LtCol values are all larger than Robinson's, 
but the estimates are all within range of one standard 
deviation and, therefore, competitive. MLE, once again, 
shows well for this rank. In opposing fashion to that of 
the combat support aggregate, the FTE inversion outperformed 
all Basic inversions for all variables. 

The MAD measure of effectiveness. Table 9 confirms 
the results found above in that Mod EB was best for IstLt's 
and MLE did well for LtCol 's. The values for IstLt were 
found to be too large overall, while all estimates for LtCol 
were competitive and in the desired range. The mix between 
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Table 8. GROUND COMBAT FIGURES OF MERIT 





1981 19S2 1983 

TRANSFORMED FOM 


% 


1st Lt 


MLE 


2.5764 


4.4136 


6.8744 




TSCA 


1.9741 


3.8072 


6.1152 




JS 


2.1737 


3.9738 


6.3970 




MOD JS 


3.0859 


4.7919 


7.5619 




EB 


2.5558 


4.3626 


6.9607 




MOD EB 


3.1538 


3.5202 


4.7069 




Lt Col 


MLE 


1.2126 


1.9674 


3.1084 




TSCA 


1.5715 


2.3553 


3.5566 




JS 


1.5571 


2.2719 


3.4365 




MOD JS 


1.6227 


2.1317 


3.1883 




EB 


1.7087 


2.2464 


3.3345 




MOD EB 


1.1964 


2.5830 

ORIGINAL FOM 


2.5029 




1st Lt 


MLE 


87.3362 


119.2637 


143.9128 




TSCA 


87.9712 


116.6767 


145.6236 




TSCAM 


83.7648 


111.3995 


142.0902 




JS 


95.4717 


122.9425 


152.7354 




JSM 


92.0238 


117.5264 


149.3965 




MOD JS 


124.6490 


147.6611 


180.5502 




MOD JSM 


134.9472 


142.5461 


178.6173 




EB 


109.6349 


136.6321 


167.3234 




EBM 


105.5708 


131.0400 


163.8189 




MOD EB 


64.3647 


81.2071 


101.9511 




MOD EBM 


63.6623 


77.8159 


100.7150 




Lt Col 


MLE 


35.0538 


45.8172 


53.0111 




TSCA 


55.9649 


198.8659 


123.1306 




TSCAM 


43.5234 


201.9953 


122.8648 




JS 


56.3346 


198.9037 


122.2786 




JSM 


43.9957 


201.3299 


115.7143 




MOD JS 


60.0730 


202.2959 


122.9959 




MOD JSM 


52.3940 


210.7551 


121.4254 




EB 


62.7441 


205.4249 


126.6544 




EBM 


50.4049 


209.3853 


125.7642 




MOD EB 


48.0871 


197.7316 


111.9375 




MOD EBM 


38.1231 


323.5654 


242.4398 





- » 
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over-estimates and under-estimates, i.e., the number of 

cells forecasted over or under the actuals, is not a 

* 

desirable one, yet the same mix holds true for all 
validation years. This also was the finding of both 

preceding aggregates. 



Table 9. MEAN ABSOLUTE DEVIATION (GROUND COMBAT) 





1981 


1982 


1983 


Isl’ 

MLE 


UNDERAGE 


0.1525( 5) 


0.1855(10) 


0.2466(10) 


OVERAGE 


1.9322(30) 


2.5979(25) 


3.1691(25) 


TOTAL 


2.0847(35) 


2.7834(35) 


3.4157(35) 


MOD EB 


UNDERAGE 


0.2340(10) 


0.1323( 8) 


0.2220( 8) 


OVERAGE 


1.7743(33) 


2.3086(32) 


2.8449(30) 


TOTAL 


2.0083(43) 


2.4409(40) 


3.0669(38) 


MOD EBM 


UNDER.AGE 


0.2618(10) 


0.1705(10) 


0.2617(10) 


OVER.AGE 


1.6492(33) 


2.1738(30) 


2.6935(28) 


TOTAL 


1.9111(43) 


2.3443(40) 


2.9552(38) 


Lt Col 
MLE 


UNDER.AGE 


0.2313(10) 


0.3034(11) 


0.1368(10) 


overage 


0.4182(25) 


0.4202(24) 


0.6453(25) 


TOTAL 


0.6495(35) 


0.7236(35) 


0.7821(35) 


MOD EB 


UNDERAGE 


0.2315(11) 


0.3299(14) 


0.1426(12) 


OVERAGE 


0.6392(34) 


0.6922(32) 


0.9523(31) 


TOTAL 


0.8706(45) 


1.0221(46) 


1.0950(43) 


MOD EBM 


UNDERAGE 


0.2650(12) 


0.3790(15) 


0.1864(12) 


OVERAGE 


0.4071(33) 


0.4427(31) 


0.6535(31) 


TOTAL 


0.6720(45) 


0.8217(46) 


0.8399(43) 
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Test Aggregate Case 



4 . 



In an effort to experiment with aggregate groupings, 
another case was built upon which all estimation methods 
were applied, and the results compared with previously 
studied groups. The grade was held fixed, in this case the 
grade of Captain was used, with a LOS in the eighth year. 
Several MOS ' s were chosen so that to insure an aggregate 



with non- 


-zero inventories. In the transformed scale, Table 


10, all 


estimates have small and competitive values, with 


Mod EB 


and MLE doing well. In the original scale, once 



again, all estimates are well within range, but Mod EB and 
MLE were the standouts. The FTE inversion option did not 
fair well with this aggregate. All estimates in the MAD 
measure of effectiveness. Table 11, are very competitive, 
and the first validation year mix (overage and underage) may 
be the best yet. It is interesting to note that the cell 



split in 


the first validation year between overage and 


underage 


forecasts is contrary to that of the historical 


displays 


above. Although the total number of cells used 


(k=6) is 


much smaller than the previous table, it is felt 


that the 


aggregate scheme was the major factor in this 


finding , 


and the results are promising. 
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Table 10, TEST AGGREGATE FIGURES OF MERIT 





1981 


1982 


1983 




TRANSFORMED FOM 




Captain 




■ 




MLE 


0.9603 


2.3463 


2.0842 


TSCA 


1.4151 


3.3033 


2.9152 


JS 


1.3753 


3.3969 


2.8725 


MOD JS 


1.3943 


3.8165 


2.8829 


EB 


1.6840 


4.2698 


3.2957 


MOD EB 


0.8249 


1.7132 


1. 2980 






ORIGINAL FOM 




Captain 








MLE 


3.0402 


4.8104 


4.9399 


TSCA 


4.6161 


10.2779 


7.5714 


TSCAM 


8.1467 


7.6241 


8.8347 


JS 


4.3970 


10.3000 


7.4726 


JSM 


7.3199 


7.6666 


8.3860 


MOD JS 


4.2990 


10.7640 


7.6371 


MOD JSM 


6.0079 


8.3070 


7.9259 


EB 


5.3660 


12.0337 


8.8679 


EBM 


6.7181 


9.3915 


8.7959 


MOD EB 


1.4463 


6.0323 


3.9066 


MOD EBM 


4.6465 


4.3354 


5.9175 



- » 
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Table 1 1. MEAN ABSOLUTE DEVIATION (TEST CASE) 



' 


1981 


1982 


1983 


Captain 








MLE 

UNDERAGE 


0.2625( 3) 


0.0463( 1) 


0.2097( 2) 


OVERAGE 


0.135l{ 3) 


0.7374( 5) 


0.4391( 4) 


TOTAL 


0.3976( 6) 


0.7837( 6) 


0.6488( 6) 


MOD EB 
UNDERAGE 


0.1342( 4) 


0.0000( 0) 


0.0965( 1) 


OVERAGE 


0.2575( 2) 


0.960 1( 6) 


0.6306{ 5) 


TOTAL 


0.39 18( 6) 


0.960 1( 6) 


0.7270( 6) 


MOD EBM 
UNDERAGE 


0.37 12( 4) 


0.0809( 1) 


0.2352( 2) 


OVEIUAGE 


0.0852( 2) 


0.6148( 5) 


0.3025( 4) 


TOTAL 


0.4563( 6) 


0.6957( 6) 


0.5377( 6) 



. * 
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V. CONCLUSIONS AND RECOMMENDATIONS 



A. RESULTS 

Omitting the fact of differences in computational 
methods between this paper and the Robinson thesis, it is 
clear that the results are still inconclusive at this time. 
In the historical comparisons, the estimates studied varied 
in their performance from one aggregate to another. The 
different risk calculations, in the transformed and original 
scales, as well as the different measures of effectiveness 
showed little agreement in estimation methods. Long range 
estimation still remains a problem since with time, the 
estimate risks increase significantly. Much emphasis was 
put on the modified empirical Bayes estimator in this paper 
and those specific results look promising under possibly 
different and more stable aggregate conditions. 

One area of mention here concerning the results has to 
do with the composition of inventories used in the selected 
aggregates. It was noticed that the range of inventory 
(smallest to largest value in any particular MOS) was quite 
large, and may be affecting the results. In fact, Carter 
and Rolph [Ref. 4:p. 883] managed this point by further 
disaggregating into levels of activity as it concerned their 
work with fire alarm estimation. 
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B. CONCLUSIONS 



The results of the previous feasibility studies and this 
paper leads to several conclusions. The chosen aggregate 
still seems to be the key element in the work of attrition 
estimation. It does seem clear that the aggregation scheme 
should be focused on fixed ranks, in narrow LOS ranges and 
multiple MOS's. This conclusion is based upon the relative 
success of the test aggregate utilized in this work, and 
upon similar success by Carter and Rolph. Additionally, the 
partitioning of inventory levels in the aggregate should 
stabilize the estimation method performance. All the 
estimation method studied to date have their strengths, but 
still, no one method of estimation excels over the others. 

C . RECOMMENDATIONS 

No estimation method studied in this paper is 
recommended for implementation at this time. In addition, 
no method is recommended for removal from consideration 
since basic aggregation problems are felt to be the largest 
contributory factor to the mixed results found here. These 
recommendations are presented for further study. 

1. The aggregation method offered by Captain Larsen 
[Ref. 3] should be examined further with irs 
relevance to the scheme in the test case used here. 

2. Utilization of the detailed data provided by NPRDC as 
soon as possible in order to expand the available 
information and open new avenues of aggregation. 

3. Re-examine the inventory mean value used in the FTE 
inverse to determine its future role in the 
estimation process. 
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APPENDIX A. 



ESTIMATION ALGORITHMS 



A. NOTATION 

The following notation is a generalization of that used 
by Robinson [Ref. 2:p. 71] as it pertains to the aggregation 
scheme used in this paper. The indexing system has been 
changed to reflect the use of officer grades to identify the 
cells . 

i = index of MOS cells in the aggregate, 
j = index of LOS cells in the aggregate, 
k = index of Grade cells in the aggregate, 
t = index of time periods in the aggregate. 
inv(i,j,k,t) = inventory with MOS i, LOS j and Grade k 
during year t. 

y(i,j,k,t) = number of attritions in cell (i,j,k) during 
year t. 

n(i,j,k,t) = maximum I y ( i , j , k , t ) , . 5 [inv ( i , j , k , t ) + 
inv (i,j,k,t + l)]t 

D = matrix identifying cells with no inventory over all 
estimation years. 

d(i,j,k) = 0, if cell is a structural zero. 
d(i,j,k) = 1, if cell is not a structural zero. 



44 



B. FREEMAM-TUKEY EXACT INVERSION ALGORITHM 



This algorithm provides a formula for the inverse of the 
Freeman-Tukey double arcsine transformation to express the 
means of the transformed values into proportions on the 
original scale. 

STEP 1. Transform the data using the Freeman-Tukey 
double arcsine transformation. 



x(i, j ,k, t) 



0.5[n(i,j,k,t)+0.5]^' 



f sini [ 



2(y(i,i,k,t) ) 
(n ( i , j , k , t ) +1 ) 



- 1 ] 



+ sin- 



. 2 (y (i, j ,k, t)+l 
Mn(i,j,k,t)+1) 



- 1 ] 1 



(A. 1} 



STEP 2 . Calculate the average cell inventory and 
transformed values over all T. 



x(i,j,k) = (1/T)Z: x(i,j,k), (A. 2) 

t 



inv(i,j,k) = (1/T) 2: inv{i,j,k). (A. 3) 

t 

STEP 3. Check for inversion values less than zero 
or greater than one. The algorithm is terminated at this 
point if these conditions are met. 

tm = X (i , j , k) / ( inv ( i , j , k) +0 . 5) 1 ^ 2 + n/2, (A. 4) 

tmi = sin- Ml/ ( inv { i , j , k) +l))i/2, {A. 5) 



If 


tm < 


tmi , 


FTE = 0.0, 


IF 


tm > 


(n-tmi ) , 


FTE = 1.0, 
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STEP 4. 



Compute the sgn function. 



If 


cos ( tm) 


> 0.0 


sgn = 1.0, 


If 


cos ( tm) 


< 0.0 


sgn = -1.0. 



STEP 5. Computer the inverted value of the 
transformed variable 



FTE 



0 . 5 ( 1-sgn 1 1- [sin ( tm) 



sin(tm) - l/sin(tm) . 
inv ( i , j , k ) 



(A. 6) 



C. TRANSFORM SCALE VARIANCE ALGORITHM 

This algorithm produces the variance of the transformed 
data within the unstable variance range where the assumption 
of normality does not hold. As discussed in Chapter II, 
this range of instability was determined to be between 1.001 
and 2.2. The method for computing the variance within this 
range is presented here. 

STEP 1. Compute the transformed scale vs;lues to be 
used for the variance. Since x(i,j,k) was computed in the 
FTE inversion algorithm of this appendix, the formulae will 
not be repeated. 

x(t) = ZI [0.5 + inv (i , j , k , t) ] 1 ‘ , (A. 7) 

iik 

= (i,j,k) = x(i,j,k) + (n/2) . (A. 3) 
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STEP 2. Compute the variance of transformed cell 
values. The formula is derived in Chapter II as a linear 
regression equation. The range of the variance is limited 
to 0.05 < var(i,j,k) ^ 1.0. 

var(i,j,k) = az(i,j,k)'’Mz{i,j,k) - 1)”-, 

1.001 < z(i,j,k) ^ 2.2 (X.9) 

where a = 1.6835, bl = -0.3934, and b2 = 0.8991. 

D. EMPIRICAL BAYES SHRINKAGE ALGORITHM 

This algorithm computes the empirical Bayes shrinkage 
values and estimator of attrition rates. Prior to entering 
this algorithm, the value for a has been set to 0.0. 

STEP 1. Initialization. 

ai=a. (A. 10) 

STEP 2. Define the a and Y parameters in the 
e.mpirical Bayes context. 

a = 21 l/(a + var(i,j,k)), (.\.ll) 

ijk 

where var(i,j,k) is defined by formula (A. 9). 

T = 21 a/ H l/(a + var(i,j,k)). (A. 12 

ijk ijk 
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STEP 3. Compute the empirical Bayes prior variance. 



k-1- Z: a[x(i,j,k)-z]2 

a = a - ijk , (A. 13) 

Z: [a (x (i , j , k) - z] = 

ijk 



where z = Z: 7 'x(i,j,k) and k = El d(i,j,k). 
ijk ijk 



STEP 4. Check for the following conditions and 
branch accordingly. 



If, a < 0.0 set a = 0.0 and go to step 5, 

If, ja - ai I > .0001 return to Step 1. 



STEP 5. Compute the empirical Bayes estimator. 



x^Mi, j,k) = [ 



a+var ( i , j , k] 



] X ( i , j , k ) + 



(A. 14) 



. var(i,j,k) 
^a+var { i , j , k) 



STEP 6. Compute the empirical Bayes shrinkage set. 

shr^®(i,j,k) = var ( i , j , k) / { a + var(i,j,k)). (A. 15) 

E. MODIFIED EMPIRICAL BAYES SHRINKAGE ALGORITHM 

This algorithm differs only slightly from that of the 
empirical Bayes shrinkage algorithm listed above, and for 
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that reason. 



only the pertinent steps will be displayed 



here . 

Consider a particular cell and its transformed value 
:■: ( t ) for t=l,...,T in the estimation set. To form a 
weighted time average, set 

xt{t) = x(t)/[0.5 + inv(t)]i ''2 (A. 15) 

This leads to a modification of step 1 in the following 
form. 

STEP 1. Compute the weighted transformed scale 

average . 

^ ^ Z xt(t) (A. 17) 

t 

Since the weights will modify the variance of the weighted 
transformed scale values. 

var[xt(t)] = var ( i , j , k ) / ( 0 . 5 + inv(t)) (A. 18) 

var[xt(t)] = (m)^ 21 var[xt(t)] = vt/T (A. 19) 

t 

The extreme right of formula A. 19 serves to define vt , 
the approximation for the weighted variance in the 
validation set for that cell. 

The values of xt and var(xt) are carried throughout this 
version of empirical Bayes estimation to form the values 
XTEB for the empirical Bayes attrition rate. For the 
purposes of inversion to the original scale, the usual Basic 
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formula is modified since the weighting factor, 0.5 + inv, 
is already incorporated into the transform value. The new 
Basic inversion for this version is simply as follows: 

p = 1/2[1 + sin(XTEB)] (A. 20) 
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APPENDIX B. MARINE CORPS PERSONNEL INVENTORY 
ATTRITION ANALYSIS PROGRAM 



PROGRAM MARCOR ^ ^ 

* MARINE CORPS PERSON’NEL INVENTORY ANALYSIS PROGRAM * 

iV * 



if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 



if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 

if 



PURPOSE: TO PROVIDE FUTURE PERSONNEL LOSS ESTIMATION OF THE * 
MARINE CORPS OFFICER STRUCTURE BY SEVERAL DIFFERENT * 
ESTIMATION SCHEMES. * 

PROGRAMMERS: LUIS URIBE, INDEPENDENT CONTRACTOR * 

CAPT C R DICKINSON, USMC * 

DESCRIPTION: WRITTEN IN FORTRAN 77 FOR THE IBM 3033 MAINFRAME * 
COMPUTER RESIDENT AT NPS, THIS PROGRAM INCORPORATES * 
METHODS EXPLORED BY MAJOR D. D. TUCKER , MAJOR J. R. * 
ROBINSON AND CAPTAIN DICKINSON. ESTIMATION SCHEMES * 



INCLUDED ARE: * 

MAXIMUM LIKELIHOOD ESTIMATION (MLE) * 

TRANSFOR.MED SCALE CELL AVERAGE (TSCA) 

JAMES -STEIN * 

LIMITED TRANSLATION JAMES-STEIN * 

EMPIRICAL BAYES ESTIMATION * 

RISK ANALYSIS WAS PERFORMED IN BOTH THE ORIGINAL AND * 
TRANSFORMED SCALES. ON THE TRANSFORMED SCALE, * 

ACHIEVED BY USING THE FREEMAN-TUKEY DOUBLE ARCSINE * 
TRANSFORMATION, THE RISK WAS DEFINED BY THE AVERAGE * 
SQUARED DEVIATION OF ACTUAL AND ESTIMATED VALUES. * 

FOR THE ORIGINAL SCALE, THE TRANSFORMED VALUES WERE * 
FIRST INVERTED USING TWO TECHNIQUES. AN AD-HOC * 

FORMULA USED BY TUCKER AND ROBINSON AND AN EXACT * 

INVERSE FORMULA FROM JOHN J MILLER. TWO MEASURES OF * 
EFFECTIVENESS, CHI-SQUARE GOODNESS OF FIT A.NT) MEAN * 
ABSOLUTE DEVIATION (MAD), WERE THEN EVALUATED. * 

RESULTS ARE PRINTED, WITH ALL ARRAYS COMPUTED, AT THE * 
COMPLETION OF THE PROGRAM. * 

VARIABLE DESCRIPTION: * 

THE FOLLOWING LIST OF VARIABLES, WITH A BRIEF * 

DESCRIPTION, FORM THE BASIS FOR THE VARIABLE NAMES * 
FOR THE ENTIRE PROGRAM. * 

MLE - MAXIMUM LIKELIHOOD ESTIMATOR * 

TSCA - TRANSFORMED SCALE CELL AVERAGE * 

J - JAMES-STEIN ESTIMATOR * 

J1 - JAMES-STEIN ESTIMATOR MODIFIED BY * 

A VARIATION ON SHRINKAGE * 

JLT - LIMITED TRANSLATION JAMES-STEIN * 

EB - EMPIRICAL BAYES * 

TEB - EMPIRICAL BAYES MODIFIED BY A * 

VARIATION ON * 

TRANSFORMED SPACE VARIABLES * 

THE FOLLOWING ARE UNIQUE VARIABLES IN THE TRANSFORMED * 
SPACE: * 

X - VALUE OF THE FREEMAN-TUKEY * 



MAROOOlO 

MAR00020 

MAR00030 

MAR00040 

MAR00050 

MAR00060 

MAR00070 

MAR00080 

MAR00090 

MAROOlOO 

MAROOllO 

MAR00120 

MAR00130 

MR00140 

MAR00150 

MAR00160 

MAR00170 

MAR00180 

MAR00190 

MAR00200 

MAR00210 

MAR00220 

MAR00230 

MAR00240 

MAR00250 

MAR00260 

MAR00270 

MAR00280 

MAR00290 

MAR00300 

MAR00310 

MAR00320 

MAR00330 

MAR00340 

MAR00350 

MAR00360 

MAR00370 

MAR00380 

MAR00390 

MAR00400 

MAR00410 

MAR00420 

MAR00430 

MAR00440 

MAR00450 

MAR00460 

MAR00470 

MAR00480 

MAR00490 
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•k 


TRANSFORMATION 


k 


MAR00500 


k 


XB - (X-BAR) AVERAGE TRANSFORM VALUE 


k 


MAR00510 


k 


OVER ESTIMATION YEARS 


k 


MAR00520 


k 


XBB - GRAND 1 


MEAN OF TRANSFORMED VALUES 


k 


MAR00530 


k 


SSE - SUM OF 


SQUARED DIFFERENCES, ERROR 


k 


' MAR00540 


k 


(X AND 


X-BAR) 


k 


MAR00550 


k 


SST - SUM OF 


SQUARED DIFFERENCES, TOTAL 


k 


MAR00560 


k 


(X AND 


XBB) 


k 


MAR00570 


k 


SSB - SUM OF 


SQUARED DIFFERENCES, BETWEEN 


k 


MAR00580 


k 




(SST-SSE) 


k 


MAR00590 


k 


SHR - 


JAMES -STEIN SHRINKAGE VALUES 


k 


MAR00600 


k 


SHRl - 


MODIFIED JAMES-STEIN SHRINKAGE 


k 


MAR00610 


k 




VALUES 


k 


MAR00620 


k 


SHREB - 


EMPIRICAL BAYES SHRINKAGE VALUES 


k 


MAR00630 


k 


SHRTEB - 


MODIFIED EMPIRICAL BAYES SHRINKAGE 


k 


MAR00640 


k 




VALUES 


k 


MAR00650 


k 


THE FOLLOWING ARE 1 


MODIFIED BASIC VARIABLES IN THE 


k 


MAR00660 


k 


TRANSFORMED SPACE. 


THIS IS DONE BY PLACING AN 'X' FOR 


k 


MAR00670 


k 


A PURE TRANSFORMED 


VARIABLE, OR AN 'R' FOR A RISK 


k 


MAR00680 


k 


VARIABLE IN FRONT ' 


OF THE BASIC VARIABLE. 


k 


MAR00690 


k 


MLE 


RML (MLE) 


k 


MAR00700 


k 


TSCA 


RTSCA 


k 


MAR00710 


k 


XJ 


RJ 


k 


MAR00720 


k ■ 


XJl 


RJl 


k 


MAR00730 


k 


XJLT 


RJLT 


k 


MAR00740 


k 


XEB 


RXEB 


k 


MAR00750 


k 


XTEB 


RTEB 


k 


MAR00760 


k 


ORIGINAL SPACE VARIABLES 




k 


MAR00770 


k 


THE FOLLOWING ARE 1 


MODIFIED BASIC VARIABLES IN THE 


k 


MAR00780 


k 


ORIGINAL SPACE. SIMILAR TO WHAT IS DONE IN THE 


k 


MAR00790 


k 


TRANSFORMED SPACE, 


A 'P' SIGNIFIES A PURE INVERTED 


k 


MAR00800 


k 


VARIABLE, AND AN 


RR' REPRESENTS A RISK VARIABLE. AN 


k 


MAR00810 


k 


'M' AT THE END OF , 


A VARIABLE NAME IDENTIFIES A MILLER 


k 


MAR00820 


k 


INVERSE VARIABLE. 




k 


MAR00830 


k 


PMLE 


RRML (MLE) 


k 


MAR00840 


k 


PJ 


RRJ 


k 


MAR00850 


k 


PJM 


RRPJM 


k 


MAR00860 


k 


PJl 


RRJl 


k 


MAR00870 


k 


PJIM 


RRJIM 


k 


MAR00880 


k 


PJLT 


RRJLT 


k 


MAR00890 


k 


PEB 


RPJflEB 


k 


MAR00900 


k 


PEBM 


RRXEBM 


k 


MAR00910 


k 


PTEB 


RRRTEB 


k 


MAR00920 


k 


PTEBM 


RRTEBM 


k 


MAR00930 


k 


THE FOLLOWING LIST 


IS OF PREFIXES USED TO DESCRIBE 


k 


MAR00940 


k 


THE MEAN ABSOLUTE ] 


DEVIATION VARIABLES: 


k 


MAR00950 


k 


MO - MEAN DEVIATION 'OVERAGE' 


k 


MAR00960 


k 


MU - MEAN DEVIATION 'UNDERAGE' 


k 


MAR00970 


k 


MD - MEAN 'TOTAL DEVIATION 


k 


MAR00980 


k 


KVO - NUMBER 


OF 'OVERAGE' ESTIMATES IN THE 


k 


MAR00990 


k 


AGGREGATE 


k 


MAROIOOO 


k 


KVU - NUMBER 


OF 'UNDERAGE' ESTIMATES IN THE 


k 


MAROlOlO 


k 


AGGREGATE 


k 


MAR01020 


k 


KVD - NUMBER 


OF 'TOTAL' ESTIMATION DEVIATIONS 


k 


MAR01030 


k 


IN THE 


AGGREGATE 


k 


MAR01040 


k 


THESE PREFIXES ARE 


APPENDED TO THE BASIC VARIABLES TO 


k 


MAR01050 



.. ^ 
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* 




FORM THE MAD VARIABLES AS FOLLOWS 






>v 


MAR01060 


if 




MOTSCA 


MUTSCA 


MDTSCA 


KVOTSA 


KVUTSA 


KVDTSA 


(TSCA) 


if 


MAR01070 


if 




MOTSAM 


MUTSAM 


MDTSAM 


KVOTSM 


KVUTSM 


KVDTSM 


(TSCAM) 


if 


MAR01080 


if 




MOPMLE 


MUPMLE 


MDPMLE 


KVOPML 


KWPML 


KVTDPML 


(PMLE) 


if 


MAR01090 


if 




MOPJ 


MUPJ 


MDPJ 


KVOPJ 


KVUPJ 


K*VDPJ 


(PJ) 


if 


' HARO 11 00 


if 




MOPJM 


MUPJM 


MDPJM 


KTOPJM 


KVUPJM 


KVDPJM 


(PJM) 


if 


MAROlllO 


if 




MOPJl 


MUPJl 


MDPJl 


KVOPJl 


KVUPJ 1 


KTDPJl 


(PJl) 


if 


MAR01120 


if 




MOPJIM 


MUPJIM 


MDPJIM 


KVOPIM 


KVUPIM 


KVDPIM 


(PJIM) 


if 


MAR01130 


if 




MOPEB 


MUPEB 


MDPEB 


KVOPEB 


KVUPEB 


KVDPEB 


(PEB) 


if 


MAR01140 


if 




MOPEBM 


MUPEBM 


MDPEBM 


KVOPBM 


K\'UPBM 


KVDPBH 


(PEBM) 


if 


MAR01150 


if 




MOPTEB 


MUPTEB 


MDPTEB 


KVOPTB 


KVUPTB 


KVDPTB 


(PTEB) 


if 


MAR01160 


if 




MOPTBM 


MUPTBM 


MDPTBM 


KVOBTM 


KVUBTM 


KVDBTM 


(PTEBM) 


if 
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if 


INPUT/OUTPUT: 














if 
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if 




INPUT IS TAKEN FROM A FILE CREATED FROM A TAPE 


if 
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if 




PROVIDED BY NPRDC 










if 
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if 




OUTPUT IS TO THE PRINTER 










if 


MAR01210 
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ifif 


FIXED 


PARAMETERS 
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PARAMETER (hD(YR=7, MXLOS=31, MXGRD=10, MXMOS=15) 






MAR01240 


ifif 


INPUT 


PARAMETERS 
















MAR01250 



INTEGER STl, ST2, LYR 
INTEGER SLOSl, SLOS2 
INTEGER SMOS(MXMOS), NMOS 
INTEGER SGRD(MXGRD), NGRD 
** VARIABLES IN ORIGINAL SPACE 

INTEGER T, FLAG, IVYR(MXYR) 

INTEGERS DCMXMOS, MXLOS, MXGRD) 

INTEGER KV(MXYR) 

REAL AVCINV,AINV,ATT,ATT1 

REAL CINV(MXMOS,MXLOS, MXGRD, MXYR) , Y(MXMOS ,MXLOS , MXGRD ,MXYR) 
REAL ACINV(MXMOS,MXLOS,hD(GRD) 

REAL P J(MXMOS,MXLOS, MXGRD ) , PJM(MXMOS ,MXLOS-, MXGRD) 

REAL PMLE(MXMOS,MXLOS, MXGRD), PJLT(MXMOS ,MXLOS , MXGRD) 

REAL PJl(MXMOS,MXLOS, MXGRD), PJ1M(MXM0S ,MXLOS , MXGRD) 

REAL PEB(MXMOS,MXLOS, MXGRD), PEBM(MXMOS ,MXLOS , MXGRD) 

REAL PTEB(MXMOS ,MXLOS , MXGRD) , PTEBM(MXMOS ,MXLOS , MXGRD) 

REAL RRTSCAC MXYR ) ,RRJ( MXYR) ,RRJLT(MXYR) ,RRML(MXYR) 

REAL RRTSAM(MXYR) ,RRJM(MXTR) ,RRJ1(MXYR) ,RRJ1M(MXYR) 

REAL RRPEB(MXYR),RRPEBM(MXYR) 

REAL RRPTEB(MXYR) ,RRPTBM(MXYR) 

** VARIABLES IN TR^NSFORMED SPACE 

REAL X(MXMOS,MXLOS, MXGRD, MXYR) 

REAL XT(MXMOS,MXLOS, MXGRD, MXYR) 

REAL XB(MXMOS,MXLOS, MXGRD), TSCA(MXMOS ,MXLOS , MXGRD) 

REAL XTB(MXMOS,MXLOS, MXGRD) 

REAL XJ(MXMOS,MXLOS, MXGRD), XJLT(MXMOS ,MXLOS , MXGRD) 

REAL XEB(MXMOS,MXLOS, MXGRD), XJl(MXMOS ,MXLOS , MXGRD) 

REAL XTEB(MXMOS,MXLOS, MXGRD) 

REAL XMLE(MXMOS,MXLOS, MXGRD) 

REAL SHREB(MXMOS,MXLOS, MXGRD) 

REAL SHRTEB(MXMOS,MXLOS, MXGRD) 

REAL RML(MXYR), RSL( MXYR ) ,RXEB( MXYR) 

REAL RXTEB(MXYR) 

REAL RTSCAC MXYR), RTJ( MXYR), RTJLTC MXYR) 

REAL RJ(MXYR) ,RJLT(MXYR) , RJl(MXYR) 

REAL TSCAM(MXMOS,MXLOS, MXGRD) 
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REAL V( MXMOS , MXLOS , MXGRD ) 

REAL VT( MXMOS, MXLOS, MXGRD) 

REAL*8 SHR, SHRX, SHRl, SSB, SST, SSE, XBB 
** VARIABLES TOR FIGURE OF MERIT 

INTEGER KVOTSA(MXYR),KVOTSM(MXYR),KVOPJ(MXYR),KVOPJl(MXYR) 

INTEGER KVOPJM(MXYR) ,KVOPEB( MXYR) ,KV0P1M(MXYR) ,KVOPBM(MXYR) 
INTEGER KVOPTB(MXTR) ,KVOBTM(MXYR) 

INTEGER KVOPML(MXYR) 

INTEGER KVUTSA(MXYR) ,KVUTSM(MXYR) ,KVUPJ(MXYR) ,KVUPJ1(MXYR) 

INTEGER KVUPJM(MXYR) ,KVUPEB(MXYR) ,KVUP1M(MXYR) ,KVUPBM(MXYR) 
INTEGER KVUPTBC MXYR ),KVUBTM( MXYR) 

INTEGER KVUPML(MXYR) 

INTEGER KVDTSA(MXYR) ,KVDTSM(MXYR) ,KVDPJ(MXYR) ,KVDPJ1(MXYR) 

INTEGER KVDPJM(MXYR) ,KVDPEB(MXYR) ,KVDP1M(MXYR) ,KVDPBM(MXYR) 
INTEGER KVDPTB(^LXYR),KVDBTM(MXYR) 

INTEGER KVDPML(MXYR) 

REAL MUTSCA(MXYR) ,MUTSAM(MXYR) ,MUPJ(MXYR) ,MUPJ1(MXYR) ,MUPJM(MXYR) 
REAL MUPJIM(MXYR) ,MUPEBM(MXYR) ,MUPEB(MXYR) ,MUPMLE(MXYR) 

REAL MUPTBM(MXY'R),MUPTEB(MXYR) 

REAL MOTSCA(^LXYR) ,MOTSAM(MXYR) ,MOPJ(MXYR) ,MOPJl(MXYR) ,MOPJM(MXYR) 
REAL MOPJIM(MXYR) ,MOPEBM( MXTR) ,MOPEB(MXYR) ,MOPMLE(MXYR) 

REAL MOPTBMC MXYR ),MOPTEB( MXYR) 

REAL MDTSCA(MXYR) ,MDTSAM(MXYR) ,MDPJ(MXYR) ,MDPJ1(MXYR) ,MDPJM(MXYR) 
REAL MDPJIM(MXYR) ,MDPEBM(ffiYR) ,MDPEB(MXYR) ,MDPMLE(MXTR) 

REAL MDPTBMC MXYR ),MDPTEB( MXYR) 

REAL MILLER 

** INPUT DATA AREAS 

INTEGER MOS, LOS, GRADE, YR, INV 
INTEGER M0S1,L0S1,GRADE1,YR1,INV1 
** INPUT INITIALIZATION 

DATA KV/MXYR*0/ 

DATA MOPMLE/MXYR*0/,MUPMLE/MXYR*0/,MDPMLE/MXYR*0/ 

DATA KVOPML/MXYR*0/ ,KVUPML/MXTR*0/ ,KVDPML/MXYR*0/ 

DATA MOTSCA/MXYR*0/ ,MUTSCA/MXYR*0/ ,MDTSCA/MXYR*0/ 

DATA KVOTSA/^LXYR*0/ , KVUTSA/MXYR*0/ , KVDTSA/^LXYR*0/ 

DATA MOTSAM/MXYR*0/ ,MUTSAM/MXYR*0/ ,MDTSAM/MXTR*0/ 

DATA KVOTSM/MXYR*0/ , KVUTSM/MXYR*0/ , KVT)TSM/MXYR*0/ 

DATA MOPJ/MXYR*0/ ,MUPJ/MXYR*0/ ,MDPJ/MXYR*0/ 

DATA KVOPJ/MXYR*0/ ,KVUPJ/MXY’R*0/ ,KVDPJ/MXYR*0/ 

DATA MOPJM/MXYR*0/ ,MUPJM/MXYR*0/ ,MDPJM/MXYR*0/ 

DATA KVOPJM/MXYR*0/ ,KVUPJM/MXYR*0/ ,KVDPJM/MXYR*0/ 

DATA M0PJ1/MXTR*0/ ,MUPJ1/MXYR*0/ ,MDPJ1/MXYR*0/ 

DATA KVOPJ1/MXYR*0/,KVUPJ1/MXYR*0/,KVDPJ1/MXYR*0/ 

DATA M0PJ1M/MXYR*0/ ,MUPJ1M/MXYR*0/ ,MDPJ1M/>KYR*0/ 

DATA KV0P1M/MXYR*0/ ,KVUP1M/MXYR*0/ ,KVDP1M/MXYR*0/ 

DATA MOPEB/MXYR*0/ ,MUPEB/MXYR*0/ ,MDPEB/MXYR*0/ 

DATA KVOPEB/MXYR*0/ ,KVUPEB/MXYR*0/ ,KVDPEB/MXYR*0/ 

DATA MOPEBM/MXYR*0/ ,MUPEBM/MXYR*0/ ,MDPEBM/MXYR*0/ 

DATA KVOPBM/MXYR*0/ ,KVUPBM/MXYR*0/ ,KVDPBM/MXYR*0/ 

DATA MOPTEB/MXYR*0/ ,MUPTEB/MXYR*0/ ,MDPTEB/MXYR*0/ 

DATA KVOPTB/MXYR*0/,KVUPTB/MXYR*0/,KVDPTB/MXYR*0/ 

DATA MOPTBM/MXYR*0/,MUPTBM/MXYR*0/ ,MDPTBM/MXYR*0/ 

DATA KVOBTM/MXYR*0/ ,KVUBTM/MXYR*0/ ,KVDBTM/MXYR*0/ 
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** READ PARAMETER DATA TO USE IN THIS RUN. ** 

WRITE(5.100) 

100 FORMATC^ ENTER 1ST YEAR AND LAST YEAR TO USE FOR ESTIMATION') 
READ(5,*) ST1.ST2 

WRITE(6,*) 'ESTIMATION YEARS: ',ST1,ST2 

* 

WRITE(5.101) 

101 FORMATC* ENTER LAST YEAR PRESENT IN THE DATA BASE’) 

READ(5,*) LYR 

WRITE(6,*) 'LAST YEAR AVAILABLE FOR VALIDATION: ’ ,LYR 

* 

WRITE(5,102) 

102 FORMATC ‘ ENTER NO. OF MOS FOLLOWED BY ARRAY OF SAME LENGTH') 
READC5,*) NMOS, (SMOS(I), I=1,NH0S) 

WRITE(6,*) 'MOS SELECTED:’, (SMOS(I), 1=1, NMOS) 

* 

WRITEC5.103) 

103 FORMATC* ENTER 1ST AND LAST LOS VALUE TO USE’) 

READC5,*) SL0S1.SL0S2 

WRITEC6,*) 'LOS RANGE; ’ ,SL0S1,SL0S2 

* 

WRITEC5.104) 

104 FORMATC * ENTER NO. OF GRADE FOLLOWED BY ARRAY OF SAME LENGTH’ ) 
READC5,*) NGRD, CSGRDCI), I=1,NGRD) 

WRITEC6,*) 'GRADES SELECTED’, CSGRDCI), 1=1, NGRD) 

* 

WRITEC5.105) 

105 FORMATC ENTER DEE FACTOR’) 

READC5,*) DEE 

WRITEC6,*) 'DEE FACTOR: ’, DEE 

* 

** COMPUTE START AND END YEARS FOR VALIDATION 

NYR=LYR-ST1+1 
NVYR=LYR-ST2 
NL0S=SL0S2-SL0S1+1 
T=ST2-ST1+1 
DO 41 I=1,NVYR 
IVYRCD=ST2+I 
41 CONTINUE 

** INITIALIZE D ARRAY TO ALLOW FINDING MISSING VALUES. ** 



DO 1 I=1,MXM0S 
DO 2 J=1,MXL0S 
DO 3 K=1,MXGRD 
DCI,J,K)=-9999 
3 CONTINUE 

2 CONTINUE 

1 CONTINUE 



** • TAPE PROCESSING TO INPUT INVENTORY AND ATTRITION VALUES, ** 

** DEVELOPS D MATRIX ACCORDING TO EXTRACTION CRITERIA. ** 

**Tfr**TVi'r*******?'r***'*******Tlr****Vc**'AT^f ■*■*******■*•*• *•*****■* 



FLAG=0 

IE0F=0 
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* 

* 

11 



■‘ir/r 

* 

12 

** 



8 



CALL READER(M0S1, LOS 1, GRADE 1,YR1,INV1,ATT1, NMOS.SMOS, 

* NGRD,SGRD, SL0S1,SL0S2, ST1,ST2, IMl, ILl , IGl, ITl , lEOF) 

IFCIEOF .ME. 0) THEN 

WRITE(6,*) ***** NO DATA MEETS SELECTIONS REQS’ 

STOP 

END IF 

IFCIEOF .NE. 0) GO TO 9 
MOS=MOSl 
LOS=LOSl 
GRADE=GRADE1 
YR=YR1 
INV=INV1 
ATT=ATT1 
IM=IM1 
IL=IL1 
IG=IG1 
IT=IT1 

IFCYR. GT. STl) FLAG = 1 

CHECK CASE WHERE DATA BEGINS PAST 1ST YR 
IFCYR .GT. STl) CINVCIM,IL,IG,IT-1)=. 5*FL0ATCINV) 

CALL READERCMOSl, LOS 1, GRADE 1,YR1,INV1,ATT1, NMOS.SMOS, 

* NGRD.SGRD, SLOS1.SLOS2, ST1,ST2, IMl, ILl, IGl, ITl, lEOF) 
ONLY FOR CHECKING PURPOSES. MARK D TO INDICATE SOME YR PRESENT 

DCIM,IL,IG)=0 



MOSl. EQ. MOS .AND. LOSl.EQ. LOS 
GRADEl.EQ. GRADE)) GO TO 8 
1 



IFC . NOT. C lEOF. EQ. 0 .AND. 

.AND. 

CENTRAL INV. FOR YRl - 
IFCYRl .GT. YR+1) THEN 

CINVC IMl , ILl , IGl , IT1-1)=. 5*FLOATC INVl) 

AINV=. 5*FL0ATC INV) 

I=YR+1 

FLAG=1 

ELSE 

AINV=. 5*FLOATC INV+INVl) 

END IF 

CINVC IM,IL,IG,IT)=AMAX1CAINV, ATT) 

YCIM,IL,IG,IT)=ATT 

YR=YR1 

INV=INV1 

ATT=ATri 

IM=IM1 

IL=IL1 

IG=IG1 

IT=IT1 

CALL READERCMOSl, L0S1,GRADE1, YRl, INVl, ATTl, NMOS,SMOS, 

NGRD,SGRD, SLOSl,SLOS2, ST1,ST2, IMl , ILl, IGl , ITl , lEOF) 
GO TO 12 
CONTINUE 

WHEN YEARS MISSING AT THE END 
IFCYR .LT. LYR) THEN 
FLAG=1 

AINV=. 5*FLOATCINV) 

ELSE 
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AINV=FLOAT(INV) 

END IF 

CINV(IM,IL,IG,IT)=ANAX1(AINV,ATT) 

Y(IM,IL,IG,IT)=ATT 

CC=CINV(IM,IL,IG,IT) 

TEMP=-1. +2. *ATT/( 1. +CC) 

TEMP1=-1. + 2.*(1.+ATT)/(1. +CC) 

CC=ABS(TEMP) 

CC1=ABS( TEMPI) 

IF(CC .GT. 1. .OR. CCl .GT. 1.) WRITE(6,*) TEMP, TEMP1=' , 

* TEMP , TEMP 1 , ' , LO S , GR , YR= ' , SMOS ( IM ) , I L , S GRD( I G ) , IT 



GO TO 11 
CONTINUE 



* CALCULATE AVERAGE CENTRAL INVENTORY AND MAXIMUM 

* LIKELIHOOD ESTIMATE. 






DO 10 IM=l,NMOS 
DO 20 IL=l,NLOS 
DO 30 IG=1,NGRD 
SUMY=0 

ACINV(IM,IL,IG)=0 
DO 35 IT=1,T 
SUMY=SUMY+Y( IM , IL , IG , IT) 

ACINV(IM,IL,IG)=ACINV(IM,IL,IG)+CINV( IM,IL,IG,IT) 

5 CONTINUE 

IF(ACINV(IM,IL,IG) . NE. 0) THEN 

PMLE( IM, IL, IG)=SUMY/ACINV( IM, IL, IG) 

ELSE 

PMLE(IM,IL,IG)=0 

ENDIF 

ACINV(IM,IL,IG)=ACINV(IM,IL,IG)/T 
DO 40 IT=1,NYR 

IF A.NY CENTRAL INVENTORY>0 (OVER EST. YRS T) THEN D=1 
IF(CINV(IM,IL,IG,IT). GT. 0 .AND. IT. LE. T) D( IM, IL, IG)=1 
IF ANY CENTRAL INVENTORY>0 (OVER VAL. YRS > T) THEN COUNT 
IF(CIN'V(IM,IL,IG,IT). GT. 0 .AND. IT. GT. T) KV( IT-T)=KV( IT-T)+1 
TRANSFORMATION OF CINV USING THE FREEMAN-TUKEY 
DOUBLE ARCS IN FORMULA. 

CALL FTT(CINV(IM,IL,IG,IT),Y(IM,IL,IG,IT),X(IM,IL,IG,IT)) 
XT(IM,IL,IG,IT)=X(IM,IL,IG,IT) / SQRT(0. 5+CINV( IM,IL, IG,IT)) 
CONTINUE 
CO.NTINUE 
CONTINUE 
CONTINUE 






CHECK FOR MISSING DATA, COMPUTE THE TRANSFORMED 
CELL AVERAGES AND GRAND MEAN. 









Iricie 



XBB=0 

KK=0 

DO 50 IM=l,NMOS 
DO 60 IL=l,NLOS 
DO 70 IG=1,NGRD 
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IF(D(IM,IL,IG) .EQ. -9999) THEN 


MAR03850 


** 


REPORT MISSING COMBINATION IM,IL,IG & CLEAR D ENTRY 


MAR03860 




D(IM,IL,IG)=0 


MAR03870 




ILD=SL0S1+IL-1 


MAR03880 




FLAG=1 


MAR03890 




END IF 


MAR03900 




KK=KK+D(IM,IL,IG) 


MAR03910 




XB(IM,IL,IG)=0 


MAR03920 




XTB(IM,IL,IG)=0 


MAR03930 




DO 80 IT=1,T 


MAR03940 




XB(IM,IL,IG)=XB(IM,IL,IG)+X(IM,IL,IG,IT) 


MAR03950 




XTB(IM,IL,IG)=XTB(IM,IL,IG)+XT(IM,IL,IG,IT) 


MAR03960 


80 


CONTINUE 


MAR03970 




XB(IM,IL,IG)=XB(IM,IL,IG)/T 


MAR03980 




XTB(IM,IL,IG)=XTB(IM,IL,IG)/T 


MAR03990 




EMPIRICAL BAYES PREPARATION ON THE TRANSFORMED SCALE. 


MAR04000 




XEB(IM,IL,IG)=XB(IM,IL,IG) 


MAR04010 




XTEB(IM,IL,IG)=XTB(IM,IL,IG) 


MAR04020 


* 




MAR04030 




XBB=XBB+D( IM, IL, IG)*XB( IM, IL, IG) 


MAR04040 


70 


CONTINUE 


MAR04050 


60 


CONTINUE 


MAR04060 


50 


CONTINUE 


MAR04070 




XBB=XBB/KK 


MAR04080 


-kif/riti 
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** 


BAYES VARIANCE ITERATION SECTION. ** 


MAR04100 




iitirkirir/rkitirkirkitititititicirlrkitkit'kin'r>'rkititirirkirkirkirkirlrkitirkirkieirk*ifirkirki^it*ititieitie 


MAR04110 




IC=0 


MAR04120 


irk 


COMPUTE TRANSFORMED SCALE CELL VARIANCE. 


MAR04130 


499 


DO 500 IM=1.NM0S 


MAR04140 




DO 510 IL^1,NL0S 


MAR04150 




DO 520 IG=1,NGRD 


MAR04160 




IF (D(IM,IL,IG).NE. 1) GOTO 520 


MAR04170 




AVSQR = 0.0 


MAR04180 




DO 522 IT=1,T 


MAR04190 


522 


AVSQR = AVSQR + SQRT(0. 5 + CINV(IM,IL,IG,IT)) 


MAR04200 




Z = XEB(IM,IL,IG) + (AVSQR/T) * 1.5708 


MAR04210 




V(IM,IL,IG) = VAR(Z) 


MAR04220 


520 


CONTINUE 


MAR04230 


510 


CONTINUE 


MAR04240 


500 


CONTINUE 


MAR04250 


* 




MAR04260 




A=0 


MAR04270 


530 


AI=A 


MAR04280 




SALPH=0 


MAR04290 




DO 531 IM=1,NM0S 


MAR04300 




DO 532 IL=1,NL0S 


MAR04310 




DO 533 IG=1,NGRD 


MAR04320 


533 


IF (D(IM,IL,IG-).EQ. 1) SALPH = SALPH + 1. /(A + V(IM,IL,IG)) 


MAR04330 


532 


CONTINUE 


MAR04340 


531 


CONTINUE 


MAR04350 


* 


< 


MAR04360 




ZB=0 


MAR04370 


irk 


COMPUTE BAYES -ALPHA AND GAMMA PARAMETERS 


MAR04380 




DO 537 IM=1,NM0S 


MAR04390 




DO 538 IL=1,NL0S 


MAR04400 






58 





DO 539 IG=1,NGRD 


MAR04410 




IF (D(IM,IL,IG).EQ. 1) THEN 


MAR04420 




ALPH = 1. /(A + V(IM,IL,IG)) 


MAR04430 




GAM = ALPH / SALPH 


MAR04440 




ZB = ZB + GAM*XB(IM,IL,IG) 


MAR04450 




ENDIF 


MAR04460 


539 


CONTINUE 


MAR04470 


538 


CONTINUE 


MAR04480 


537 


CONTINUE 


MAR04490 




F=0 


MAR04500 




G^O 


MAR04510 




DO 540 IM=l,NMOS 


MAR04520 




DO 541 IL=l,NLOS 


MAR04530 




DO 542 IG=1,NGRD 


MAR04540 




IF (D(IM,IL,IG).EQ. 1) THEN 


MAR04550 




ALPH = 1. /(A + V(IM,IL,IG)) 


MAR04560 




F = F + (ALPH * (XB(IM,IL,IG) - ZB)**2) 


MAR04570 




G = G + (ALPH * (XB(IM,IL,IG)-ZB))**2 


MAR04580 




ENDIF 


MAR04590 


542 


CONTINUE 


MAR04600 


541 


CONTINUE 


MAR04610 


540 


CONTINUE 


MAR04620 


■Mt 


COMPUTE THE PRIOR VARIANCE (A) 


MAR04630 




A = A - (KK-1-F)/G 


MAR04640 


* 




MAR04650 




IF (A. LE. 0. ) THEN 


MAR04660 




A=0 


MAR04670 




GOTO 550 


MAR04680 




ENDIF 


MAR04690 




IF (ABS(A-AI).GT. . 0001) GOTO 530 


MAR04700 




COMPUTE THE EMPIRICAL BAYES ATTRITION RATES 


MAR04710 


550 


DO 570 IM=l,NMOS 


MAR04720 




DO 580 IL=l,NLOS 


MAR04730 




DO 590 IG=1,NGRD 


HAR04740 




IF (D(IM,IL,IG).EQ. 1) THEN 


MAR04750 




DEN = A + V(IM,IL,IG) 


MAR04760 




XEB(IM,IL,IG) = (A/DEN)*XB(IM,IL,IG)+((DEN“A)/DEN )*ZB 


MAR04770 




ENDIF 


MAR04780 


590 


CONTINUE 


MAR04790 


580 


CONTINUE 


MAR04800 


570 


CONTINTJE 


MAR04810 




IC=IC+1 


MAR04820 




IF (A, EQ. 0) GOTO 599 


MAR04830 




IF (IC. LT. 10) GOTO 499 


MAR04840 


* 




MAR04850 


599 


CONTINUE 


MAR04860 


* 




MAR04870 


** 


COMPUTE THE EMPIRICAL BAYES SHRINKAGE RATES 


MAR04880 




DO 600 IM=l,NMOS 


MAR04890 




DO 601 IL=l,NLOS 


MAR04900 




DO 602 IG=1,NGRD 


MAR04910 




SHREB(IM,IL,IG) = 0 


MAR04920 




IF (D(IM,IL,IG).EQ. 1) THEN 


MAR04930 




SHREB(IM,TL,IG) = V( IM, IL, IG)/(A + V(IM,IL,IG)) 


MAR04940 




ENDIF 


MAR04950 


602 


CONTINUE 


MAR04960 



59 



601 

600 

** 

** 

699 



725 



720 

710 

700 

* 

730 



733 

732 

731 



739 

738 

737 



CONTINUE 
CONTINUE 

ALTERNATE BAYES VARIANCE ITERATION SECTION. ** 

COMPUTE TRANSFORMED SCALE CELL VARIANCE. 

DO 700 IM=1,NM0S 
DO 710 IL=1,NL0S 
DO 720 IG=1,NGRD 
VT(IM,IL,IG) = 0. 0 
IF (D(IM,IL,IG).NE. 1) GOTO 720 
DO 725 IT=1,T 

ZT = X(IM,IL,IG,IT) + 1.5708 * SQRT(0.5 + 
CINV(IM,IL,IG,IT)) 

VT(IM,IL,IG) = VT(IM,IL,IG) + VAR(ZT)/(0.5 + 
CINV(IM,IL,IG,IT)) 

VT(IM,IL,IG) = VT(IM,IL,IG)/T**2 
CONTINUE 
CONTINUE 
CONTINUE 

AT=0 

ATI=AT 

SALPHT=0 

DO 731 IM=1,NM0S 
DO 732 IL=1,NL0S 
DO 733 IG=1,NGRD 

IF (D(IM,IL,IG).EQ. 1) SALPHT= SALPHT+ 1. /(AT+ VT(IM,IL,IG)) 
CONTINUE 
CONTINUE 
ZTB=0 

COMPUTE BAYES ALPHA AND GAMMA PARAMETERS 
DO 737 IM=1,NM0S 
DO 738 IL=1,NL0S 
DO 739 IG=1,NGRD 

IF (D(IM,IL,IG). EQ. 1) THEN 
ALPHT= 1. /(AT+ VT(IM,IL,IG)) 

GAMT = ALPHT / SALPHT 

ZTB = ZTB+ GAMT*XTB(IM,IL,IG) 

ENDIF 

CONTINUE 

CONTINUE 

CONTINUE 

PZTB=0. 5*( 1. +SIN( ZTB) ) 

FT=0 

GT=0 

DO 740 IM=1,NM0S 
DO 741 IL=1,NL0S 
DO 742 IG=1,NGRD 

IF (D(IM,IL,IG).EQ. 1) THEN 
ALPHT= l./(AT + VT(IM,IL,IG)) 



742 

741 



FT 

GT 

ENDIF 

CONTINUE 

CONTINUE 



= FT + (ALPHT * (XTB(IM,IL,IG) - ZTB)**2) 
= GT + (ALPHT * (XTB(IM,IL,IG) - ZTB))**2 
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MAR04980 
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MAR05000 
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MAR05020 

MAR05030 

MAR05040 
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MAR05070 

MAR05080 

MAR05090 
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MAR05130 

MAR05140 

MAR05150 

MAR05160 

MAR05170 

MAR05180 

MAR05190 

MAR05200 

MAR05210 

MAR05220 

MAR05230 

MAR05240 

MAR05250 

MAR05260 

MAR05270 

MAR05280 

MAR05290 

MAR05300 

MAR05310 

MAR05320 

MAR05330 

MAR05340 

MAR05350 

MAR05360 

MAR05370 

MAR05380 

MAR05390 

MAR05400 

MAR05410 

MAR05420 

MAR05430 

MAR05440 

MAR05450 

MAR05460 

MAR05470 

MAR05480 

MAR05490 

MAR05500 

MAR05510 

MAR05520 
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740 COxNTINUE 

** COMPUTE THE PRIOR VARIANCE (A) 

AT= AT- (KK-1-FT)/GT 

* 

IF (AT. LE. 0. ) THEN 
AT=0 

GOTO 750 
ENDIF 

IF (ABS(AT-ATI).GT. .0001) GOTO 730 
** COMPUTE THE EMPIRICAL BAYES ATTRITION RATES 

750 DO 770 IM=l,NMOS 

DO 780 IL=l,NLOS 
DO 790 IG=1,NGRD 

IF (D(IM,IL,IG).EQ. 1) THEN 
DENT = AT + VT(IM,IL,IG) 

XTEB(IM,IL,IG) = (AT/DENT)*XTB( IM.il, IG) + 

* ((DENT-AT)/DENT)*ZTB 

ENDIF 

790 CONTINUE 

780 CONTINUE 

770 CONTINUE 

** COMPUTE THE EMPIRICAL BAYES SHRINKAGE RATES 

DO 800 IM=1,NM0S 
DO 801 IL=1,NL0S 
DO 802 IG=1,NGRD 

SHRTEB(IM,IL,IG) = 0 
IF (D(IM,IL,IG).EQ. 1) THEN 

SHRTEB(IM,IL,IG) = VT( IM, IL, IG)/( AT + VT( IM, IL, IG) ) 

ENDIF 

802 CONTINUE 

801 CONTINUE 

800 CONTINUE 

** JAMES -STEIN SHRINKAGE RATES. ** 

SST=0 

SSE=0 

SSB=0 

** COMPUTE THE SUMS OF SQUARES (SST, SSE AND SSB) 

DO 90 IM=1,NM0S 
DO. 95 IL=1,NL0S 
DO 110 IG=1,NGRD 
TEMP1=0 
TEMP2=0 
DO 120 IT=1,T 

TEMP1=TEMP1+(X( IM , IL, IG, IT) -XBB)**2 
TEMP2=TEMP2+(X( IM , IL, IG, IT) -XB( IM, IL, IG) )**2 
120 CONTINUE 

SST=SST+D(IM,IL,IG)*TEMP1 
SSE=SSE+D( IM, IL, IG)*TEMP2 
110 CONTINUE 

95 CONTINUE 
90 CONTINUE 

SSB=SST-SSE 

** COMPUTE THE JAMES-STEIN SHRINKAGE RATES 

SHRXKSSE/SSB)* (KK-3)/(KK*(T-l)+2) 
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SHR=DMIN1(SHRX, 1. DO) 
SHR=DMAX1(SHR , O.DO) 
SHR1=T’'^SHR 



** LIMITED TRANSLATION OPTION & JAMES-STEIN ** 

** ON TRANSFORMED SCALE. ** 



NN=KK*(T-1) 

DO 125 IM=l,NMOS 
DO 130 IL=1,NL0S 
DO 140 IG=1,NGRD 
XJ( IM,IL,IG)=0 
XJ1(IM,IL,IG)=0 
IF (D(IM,IL,IG).EQ. 1) THEN 
XJ( IM, IL, IG)=XBB+( 1-SHR)*(XB( IM, IL, IG) -XBB) 
XJ1(IM,IL,IG)=XBB+(1-SHR1)*(XB(IM,IL,IG)-XBB) 

XJLT( IM , IL , IG)=XBB+( 1 -SHR*RHO( XB( IM,IL,IG),XBB,SSE, 

SSB,DEE,KK,T,NN))*(XB(IM,IL,IG)-XBB) 
XMLE(IM,IL,IG)=SQRT(.5+ACINV(IM,IL,IG)) 

*ASIN(2*PMLE(IM,IL,IG)-1) 

ENDIF 
140 CONTINUE 
130 CONTINUE 
125 CONTINUE 

** VALIDATION. RISK E VALUAT ED IN TRANSFORM SPACE. ** 

DO 150 IT=1, NVYR 
RTSCA(IT)=0 
RJ(IT)=0 
RJ1(IT)=0 
RXEB(IT)=0 
RXTEB(IT)=0 
C *** RJLT(IT)=0 
RML( IT)=0 
DO 160 IM=l,NMOS 
DO 170 IL=l,NLOS 
DO 180 IG=1,NGRD 

IF(D(IM,IL,IG).EQ. 1 .AND. CINV( IM, IL, IG, IT+T) . NE. 0) THEN 
RTSCAC IT)=RTSCA( IT)+( X( IM, IL, IG,T+IT) - XB( IM , IL , IG) )**2 
RJ(IT)= RJ(IT)+(X(IM,IL,IG,T+IT)- XJ( IM, IL, IG) )**2 
RJ1(IT)= RJ1(IT)+(X(IM,IL,IG,T+IT)- XJ1( IM, IL, IG) )**2 
RXEB( IT)=RXEB( IT)+(X( IM, IL, IG,T+IT) - XEB( IM, IL, IG) )**2 
RXTEB( IT)=RXTEB( IT)+( ( X( IM , IL, IG ,T+IT) /SQRT( 0.5+ 

* CINV(IM,IL,IG,IT+T))- XTEB( IM, IL, IG) )**2) 

* /VT(IM,IL,IG) 

C *** RJLT (IT)= RJLT(IT)+(X(IM,IL,IG,T+IT)-XJLT(IM,IL,IG))**2 

RML (IT)= RML ( IT)+(X( IM, IL, IG,T+IT) -XMLE( IM, IL, IG) )+*2 
ENDIF 

180 CONTINUE 

170 CONTINUE 
160 CONTINUE 

RTSCAC IT)=RTSCA( IT) /KV( IT) 

RJ (IT)=RJ (IT)/KV(IT) 

RJl (IT)=RJ1 (IT)/KV(IT) 

RXEBC IT)=RXEB( IT) /KV( IT) 
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MAR06110 
MAR06120 
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MAR06140 
MAR06150 
MAR06160 
MAR06170 
MAR06180 
MAR06190 
MAR06200 
MAR06210 
MAR06220 
MAR06230 
MAR06240 
MAR06250 
MAR06260 
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MAR06280 
MAR06290 
MAR06300 
MAR06310 
MAR06320 
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MAR06340 
MAR06350 
MAR06360 
MAR06370 
MAR06380 
MAR06390 
MAR06400 
MAR06410 
MAR06420 
MAR06430 
MAR06440 
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MAR06590 
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RXTEB( IT)=RXTEB( IT)/(KV( IT)*T) 

C *** RJLT(IT)=RJLT(IT)/KV(IT) 

RML (IT)=RML (IT)/KV(IT) 

RTJ( IT)=RTSCA( IT)/RJ( IT) 

C *** RTJLT(IT)=RTSCA(IT)/RJLT(IT) 

C *** RSL(IT)=(RJ(IT)-RTSCA(IT))/(RJLT(IT)-RTSCA(IT)) 

150 CONTINUE 

** DO RISK ANALYSIS ON THE ORIGINAL SCALE. ** 

»V aV^ V T<r*^ V iHVVr^VVcyr*iV^'iV'^^V'»V':fr*^V':frVr'fr^^VT»r‘frTV'iV^Vyryr'iV'iWrTV^V'iVVr'iVVr^V'>Wr‘iV'A'yf^V*^V^VVr'iV^Wr'A'Vr')VriV'iV'ilryr‘iV'A'yf'iV 



DO 200 IM=l,NMOS 
DO 210 IL=l,NLOS 
DO 220 IG=1,NGRD 

** AVERAGE CENTRAL INVENTORY OVER ESTIMATION SET 

CINVS=0 
DO 230 IT=1,T 

AVCINV=AVCINV + CINV(IM,IL,IG,IT) 

230 CONTINUE 

AVCINV=AVCINV/T 

TSCA(IM,IL,IG)=D(IM,IL,IG)*SCINV( XB( IM, IL, IG) , AVCINV) 

PJ(IM,IL,IG)=D(IM,IL,IG)*SCINV( XJ( IM, IL, IG) , AVCINV) 

PJ1(IM,IL,IG)=D(IM,IL,IG)*SCINV( XJ1( IM, IL, IG) , AVCINV) 

PEB(IM,IL,IG)=D(IM,IL,IG)*SCINV( XEB( IM, IL, IG) , AVCINV) 

PTEB(IM,IL,IG)=D(IM,IL,IG)*0. 5*(1. + SIN(XTEB( IM, IL, IG) ) ) 
C *** PJLT(IM,IL,IG)=D(IM,IL,IG)*SCINV(XJLT(IM,IL,IG), AVCINV) 

IF (D(IM,IL,IG).EQ. 1) THEN 

TSCAM(IM,IL,IG)=MILLER( XB(IM,IL,IG) , AVCINV ) 
PJM(IM,IL,IG)=MILLER( XJ( IM, IL, IG) , AVCINV ) 
PJ1M(IM,IL,IG)=MILLER( XJ1( IM , IL, IG) , AVCINV ) 
PEBM(IM,IL,IG)=MILLER( XEB( IM , IL, IG) , AVCINV ) 

XMIL = XTEB(IM,IL,IG) * SQRT(0. 5 + AVCINV) 
PTEBM(IM,IL,IG)=MILLER( XMIL , AVCINV ) 

END IF 

220 CONTINUE 
210 CONTINUE 
200 CONTINUE 



Vr* 



VALIDATION. RISK EVALUATED IN ORIGINAL SPACE. 
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DO 240 IT=1, NVYR 
RRTSCAC IT)=0 
RRTSAM(IT)=0 
RRJ( IT)=0 
RRJM( IT)=0 
RRPEB(IT)=0 
RRPTEB(IT)=0 
RRPEBM(IT)=0 
RRPTBMC IT)=0 
RRJ1(IT)=0 
RRJ1M(IT)=0 
C *** RRJLT(IT)=0 
RRML( IT)=0 
KP=0 

DO 250 IM=1,NM0S 
DO 260 IL=1,NL0S 
DO 270 IG=1,NGRD 
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MAR06660 
MAR06670 
MAR06680 
' MAR06690 
MAR06700 
MAR06710 
MAR06720 
MAR06730 
MAR06740 
MAR06750 
MAR06760 
MAR06770 
MAR06780 
MAR06790 
MAR06800 
MAR06810 
MAR06820 
MAR06830 
MAR06840 
MAR06850 
MAR06860 
MAR06870 
MAR06880 
MAR06890 
MAR06900 
MAR06910 
MAR06920 
MAR06930 
MAR06940 
MAR06950 
MAR06960- 
MAR06970 
MAR06980 
MAR06990 
MAR07000 
MAR07010 
MAR07020 
MAR07030 
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MAR07050 
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MAR07070 
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MAR07100 
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MAR07140 
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* 



* 



** 



* 

* 

★ 

* 

* 

* 

* 

* 

* 

* 



c *** 
c * 
** 



* 

* 



* 



* 



* 



•* 



CC=CINV( IM, IL, IG.T+IT) 

IF(CC .NE. 0) THEN 

AR=Y(IM,IL,IG,T+IT)/CC 



ELSE 



AR=0 

END IF 

IF(D(IM,IL,IG) .NE. 0 .AND. PMLE( IM, IL,IG). NE. 0 
.AND. PMLE(IM,IL,IG).NE. 1) THEN 

KP=KP+1 

RRML( IT)=RRML( IT)+( CC*( AR-PMLEC IM, IL, IG) )**2) 

/(PMLE(IM,IL,IG)*(1. - PMLE(IM,IL,IG))) 

ENDIF 

PRODUCE CHI-SQUARE MEASURE OF EFFECTIVENESS 

RRTSCAC IT)= RRTSCAC IT)+CHISQR(D( IM, IL, IG) ,TSCA( IM, IL, IG) , 

CC,AR) 

RRTSAMC IT)= RRTSAMC IT)+CHISQR( D( IM , IL, IG) ,TSCAM( IM, IL, IG) , 

CC,AR) 

RRJ (IT)= RRJ (IT)+CHISQR(D(IM,IL,IG),PJ (IM,IL,IG), 

CC,AR) 

RRJM (IT)= RRJM (IT)+CHISQR(D(IM,IL,IG) ,PJM (IM,IL,IG), 

CC,AR) 

RRJl (IT)= RRJl (IT)+CHISQR(D(IM,IL,IG),PJ1 (IM,IL,IG), 

CC,AR) 

RRJIM (IT)= RRJIM ( IT)+CHISQR(D( IM, IL, IG) ,PJ1M( IM, IL, IG) , 

CC,AR) 

RRPEB (IT)= RRPEB ( IT)+CHISQR( D( IM , IL, IG) ,PEB (IM,IL,IG), 

CC,AR) 

RRPTEB (IT)= RRPTEB ( IT)+CHISQR(D( IM, IL, IG) ,PTEB (IM,IL,IG), 

CC,AR) 

RRPEBMC IT)= RRPEBMC IT)+CHISQR(D( IM , IL, IG) ,PEBM( IM, IL, IG) , 

CC,AR) 

RRPTBMC IT)= RRPTBM(IT)+CHISQR(D( IM, IL, IG) ,PTEBM( IM, IL, IG) , 

CC,AR) 

RRJLT (IT)= RRJLT (IT)+CHISQR(D( IM, IL,IG) ,PJLT( IM, IL, IG) , 

CC,AR) 

PRODUCE MEAN ABSOLUTE DEVIATION MEASURE OF EFFECTIVENESS 

IF(D(IM,IL,IG) .NE. 0 .AND. PMLE( IM, IL, IG) . NE. 0 
.AND. PMLE(IM,IL,IG).NT:. 1) THEN 

CALL UOM( CC, AR, PMLE( IM, IL, IG) , MOPMLE(IT), 

MUPMLE(IT), KVOPML(IT), KVUPML(IT), KVDPML(IT)) 

ENDIF 



IF(CC.NE. 0. ) THEN 

CALL UOM( CC, AR, TSCA( IM,IL, IG) , MOTSCA(IT), 

MUTSCA(IT), KVOTSA(IT), KVUTSA(IT), KVDTSA(IT)) 

CALL UOM( CC, AR, TSCAM( IM, IL, IG) , MOTSAM( IT) , 

MUTSAM(IT), KVOTSM(IT), KVUTSM(IT), KVDTSM(IT)) 

CALL UOM( CC, AR, PJ( IM, IL, IG) , MOPJ(IT), 

MUPJ(IT), KVOPJ(IT), KVUPJ(IT), KVDPJ(IT)) 

CALL UOM( CC, AR, PJM( IM, IL, IG) , MOPJM(IT), 

• MUPJM(IT), KVOPJM(IT), KVUPJM(IT), KVDPJM(IT)) 

CALL UOM( CC, AR, PJ1( IM, IL, IG) , MOPJl(IT), 



MAR07210 
MAR07220 
MAR07230 
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MAR07270 
MAR07280 
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MAR07300 
MAR07310 
MAR07320 
MAR07330 
MAR07340 
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MAR07390 
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MAR07410 
MAR07420 
MAR07430 
MAR07440 
MAR07450 
MAR07460 
MAR07470 
MAR07480 
MAR07490 
MAR07500 
MAR07510 
MAR07520 
MAR07530 
MAR07540 
MAR07550 
MAR07560 
MAR07570 
MAR07580 
MAR07S90 
MAR07600 
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270 

260 

250 



240 

** 



MUPJl(IT), KVOPJl(IT), KVUPJl(IT), KVDPJl(IT)) 

CALL UOM( CC, AR, PJ1M( IM, IL, IG) , MOPJIM(IT), 

MUPJIM(IT), KVOPIM(IT), KWJPIM(IT), ICv’DriMC IT) ) 

CALL UOM( CC, AR, PEB( IM , IL, IG) , MOPEB(IT), 

MUPEB(IT), KVOPEB(IT), KVUPEB(IT), KVDPEB(IT)) 



* 



CALL UOM( CC, AR, PTEB( IM , IL, IG) , MOPTEB(IT), 

MUPTEB(IT), KVOPTB(IT), KVUPTB(IT), KVDPTB(IT)) 



CALL UOM( CC, AR, PEBM( IM, IL, IG) , MOPEBM( IT) , 

* MUPEBM(IT), KVOPBM(IT), KVUPBM(IT), KVDPBM(IT)) 

CALL UOM( CC, AR, PTEBM( IM, IL, IG) , MOPTBM(IT), 

* MUPTBM(IT), KVOBTM(IT), KVUBTM(IT), KVDBTM(IT)) 
ENDIF 

CONTINUE 

CONTINUE 

CONTINUE 

IF (KP .EQ. 0) THEN 
FF = 0 

ELSE 

FF=FLOAT( KV( IT) ) /FLOAT( KP ) 

ENDIF 

RRML( IT)=RRML( IT)*FF 
CONTINUE 

COMPUTE TIME AVERAGE MAD 

MEASURE OF EFFECTIVENESS 
DO 300 IT=1, NVYR 

MUPMLE(IT) = MUPMLE(IT)/KV(IT) 

MUTSCA(IT) = MUTSCA(IT)/KV(IT) 

MUTSAM(IT) = MUTSAM(IT)/KV(IT) 

MUPJ (IT) = MUPJ (IT)/KV(IT) 

MUPJM (IT) = MUPJM(IT)/KV(IT) 

MUPJl (IT) = MUPJ1(IT)/KV(IT) 

MUPJIM(IT) = MUPJ1M(IT)/KV(IT) 

MUPEB (IT) = MUPEB(IT)/KV(IT) 

MUPTEB (IT) = MUPTEB(IT)/KV(IT) 

MUPEBM(IT) = MUPEBM(IT)/KV(IT) 

MUPTBM(IT) = MUPTBM(IT)/K'V(IT) 

MOPMLE(IT) = MOPMLE(IT)/KV(IT) 

MOTSCA(IT) = MOTSCA(IT)/KV(IT) 

MOTSAM(IT) = MOTSAM(IT)/KV(IT) 

MOPJ (IT) = MOPJ (IT)/KV(IT) 

MOPJM (IT) = MOPJM(IT)/KV(IT) 

MOPJl (IT) = M0PJ1(IT)/KV(IT) 

MOPJIM(IT) = M0PJ1M(IT)/KV(IT) 

MOPEB (IT) = MOPEB(IT)/KV(IT) 

MOPTEB (IT) = MOPTEB(IT)/KV(IT) 

MOPEBM(IT) = MOPEBM(IT)/KV(IT) 

MOPTBM(IT) = MOPTBM(IT)/KV(IT) 

MDPMLE(IT) = MUPMLE(IT) + MOPMLE( IT) 

MDTSCA(IT) = MUTSCA(IT) + MOTSCA( IT) 

MDTSAM(IT) = MUTSAM(IT) + MOTSAM(IT) 
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300 



MDPJ (IT) = 
MDPJM (IT) = 
MDPJl (IT) = 
MDPJIM(IT) = 
MDPEB (IT) = 
MDPTEB (IT) = 
MDPEBM(IT) = 
MDPTBM(IT) = 
CONTINUE 



MUPJ ( IT) + 
MUPJM (IT) + 
MUPJl (IT) + 
MUPJIM(IT) + 
MUPEB (IT) + 
^ MUPTEB (IT) 
MUPEBM(IT) + 
NUPTBM(IT) + 



MOPJ (IT) 
MOPJM (IT) 
MOPJl (IT) 
MOPJIM(IT) 
MOPEB (IT) 

+ MOPTEB (IT) 
MOPEBM( IT) 
MOPTBM( IT) 






PRINT SECTION FOR ALL OUTPUT. 









905 



** 



906 



907 



PRINT GROUP CHARACTERISTICS 

WRITE(6,'(1X,"TEN OR LESS BAYES VARIANCE ITERATIONS")') 
WRITE(6,'(1X,I2," BAYES VARIANCE ITERATION(S) " ) ’ ) IC 
WRITE(6,'(1X,"XBB ",F12. 6)') XBB 
WRITE(6,905) T,KK,SSB,SSE,SHRX,SHR,ZB,ZTB,PZTB,A.AT 
FORMAT(' GROUP CHARACTERISTICS*, /,’ T=',I6, /,* K=',I6, 

* /,' SSB=’,F12. 4, /,' SSE=’,F12. 4, /,' SHRX=' ,F12. 4, 






/,' SHR=’,F12. 4, /,' ZB=’,F12. 4 

/,' PZTB=’ ,F12.4, /,' A=’,F12.4, /, 

PRINT RESULTS FOR TRANSFORM SPACE 
WRITE(6,*) ’ RESULTS FOR TRANSFORM SPACE’ 
WRITE(6,906) (IVYR(I), 1=1, NVYR) 

F0RMAT(14X, 10112) 

VRITE(6,907) 'RTSCA ', (RTSCA(IT), IT=1,N^/YR) 
FORMAT(* ’ ,A14,10F12. 4) 



ZTB=' ,F12.4, 
AT=’ ,F12.4) 



C 
C 
C 

•Mf 



WRITE(6,907) 'RJ 
WRITE(6,907) 'RJl 
WRITE(6,907) 'RXEB 
WRITE(6,907) 'RXTEB 
*** WRITE(6,907) 'RJLT 
WRITE( 6,907) 'RML 
WRITE(6,907) 'RSL 
*** WRITE(6,907) 'RTJ 
*** WRITE(6,907) 



(RJ (IT), 
(RJl (IT), 
(RXEB (IT), 



IT=1,NVYR) 

IT=1,NVYR) 

IT=1,NVYR) 



(RXTEB (IT), IT=1,NVYR) 
(RJLT (IT), IT=1,NVYR) 



(RML (IT), 
(RSL (IT), 
(RTJ (IT), 
(RTJLT(IT), 



’ rXJLT 

PRINT RESULTS FOR ORIGINAL SPACE 
WRITE(6,*) ' ' 

RESULTS FOR ORIGINAL SPACE’ 



IT=1,NVYR) 

IT=1,NVYR) 

IT=1,NVYR) 

IT=1,NVYR) 



WRITE(6,*) 

WRITE(6,907) ’RRML 
VRITE(6,907) ’RRTSCA 
WRITE(6,907) ’RRTSAM 
WRITE(6,907) ’RRJ 
WRITE(6,907) ’RRJM 
WRITE(6,907) ’RRJl 
WRITE(6,907) ’RRJIM 
WRITE(6,907) ’RRPEB 
WRITE(6,907) ’RRPTEB 
WRITE(6,907) ’RRPEBM 
WRITE(6,907) ’RRPTBM 
C *** WRITE(6,907) ’RRJLT 
WRITE(6,908) ’KV 
908 FORMATC ’,A14, 10112) 
VRITE(6,’(9X,’’KP 
WRITE(6,*) ' ’ 

** PRINT MAD CALCULATIONS 



(RRML (IT), IT=1,NVYR) 
(RRTSCA(IT), IT=1,NVYR) 
(RRTSAM(IT), IT=1,NVYR) 
(RRJ (IT), IT=1,NVYR) 
(RRJM (IT), IT=1,NVYR) 
(RRJl (IT), IT=1,NVYR) 
(RRJIM (IT), IT=1,NVYR) 
(RRPEB(IT), IT=1,NVYR) 
(RRPTEB(IT), IT=1,NVYR) 
(RRPEBM(IT), IT=1,NVYR) 
(RRPTBM(IT), IT=1,NVYR) 
(RRJLT(IT), IT=1,NVYR) 
(KV (IT), IT=1,NVYR) 

,112)’) KP 
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MAR08330 

MAR08340 

MAR08350 

MAR08360 

MAR08370 

MAR08380 

MAR08390 

MAR08400 

MAR08410 

MAR08420 

MAR08430 

MAR08440 

MAR08450 

MAR08460 

MAR08470 

MAR08480 

MAR08490 

MAR08500 

MAR08510 

MAR08520 

MAR08530 

MAR08540 

MAR08550 

MAR08560 

MAR08570 

MAR08580 

MAR08590 

MAR08600 

MAR08610 

MAR08620 

MAR08630 

MAR08640 

MAR08650 

MAR08660 

MAR08670 

MAR08680 

MAR08690 

MAR08700 

MAR08710 

MAR08720 

MAR08730 

MAR08740 

MAR08750 

MAR08760 

MAR08770 

MAR08780 

MAR08790 

MAR08800 

MAR08810 

MAR08820 

MAR08830 

MAR08840 

MAR08850 

MAR08860 

MAR08870 

MAR08880 



** 



1005 

1006 
1007 

** 



1015 

1016 

** 



WRITE( 6 , ' ( 2X , " PMLE " ) ' ) 

CALL WRMAD( I VYR , NVYR , MUPMLE , KVUPML , MOPMLE ,KVOPML , MDPMLE , KVDPML) 
WRITE(6,'(2X,"TSCA' ')') 

CALL WRMADC I VYR , NVYR . MUTSCA , KVUTSA ,MOTSCA , KVOTSA ,MDTSCA , KVDTSA) 
WRITEC 6 , ' ( 2X, " TSCAM ' ' ) ' ) 

CALL WRMADC I VYR , NVYR MUTS AM , KVUTSM , MOTS AM , KVOTSM , MDTSAM , KVDTSM ) 
WRITE(6,'(2X,' 'PJ' ')') 

CALL WRMADC I VYR , NVYR , MUP J , KVUP J , MOP J , KVOP J , MDP J , KVDP J ) 

WRITEC 6, 'C2X, "PJM' ')') 

CALL WRMADC I VYR , NVYR , MUP JM , KVUP JM , MOP JM , KVOP JM , MDP JM , KVDP JM ) 
WRITEC6,'C2X,"PJ1")') 

CALL WRMADC I VYR , NVYR , MUP J 1 , KVUP J 1 , MOP J 1 , KVOP J 1 , MDP J 1 , KVDP J 1 ) 
WRITEC6,'C2X,"PJ1M' ^)') 

CALL WRMAD C I VYR , NVYR , MUP J IM , KVUP IM , MOP J IM , KVOP IM , MDP J IM , KVDP IM ) 
WRITEC6,’C2X,' ’PEB")') 

CALL WRMADC I VYR , NVYR , MUPEB , KVUPEB , MOPEB , KVOPEB , MDPEB , KVDPEB ) 
WRITEC6,'C2X,"PEBM' ')’) 

CALL WRMADC I VYR , NVYR , MUPEBM , KVUPBM , MOPEBM , KVOPBM , MDPEBM , KVDPBM) 
WRITEC 6, 'C2X, "PTEB' ')' ) 

CALL WRMADC I VYR , NVYR , MUPTEB , KVUPTB , MOPTEB , KVOPTB , MDPTEB , KVDPTB ) 
WRITEC6,'C2X,"PTEBM")') 

CALL WRMADC I VYR , NVYR , MUPTBM , KVUBTM , MOPTBM , KVOBTM , MDPTBM , KVDBTM) 
WRITEC 6,*) ' ' 

PRINT THE CENTRAL INVENTORY MATRIX 
WRITEC6,*) 'CINV MATRIX' 

DO 1006 IT=1,NYR 
WRITEC6,*) ‘year’, IT 
DO 1005 IL=l,NLOS 

WRITEC 6, 1007) CCINVCIM,IL,1, IT), IM=l,NMOS) 

CONTINUE 

CONTINUE 

FORMATC 15C1X.F7. 2)) 

WRITEC 6,*) ' ' 

PRINT THE ATTRITION MATRIX 
WRITEC 6,*) 'Y MATRIX' 

DO 1016 IT=1,NYR 
WRITEC6,*) 'YEAR', IT 
DO 1015 IL=l,NLOS 

WRITEC 6, 1007) CYCIM,IL,1, IT), IM=l,NMOS) 

CONTINUE 
CONTINUE 
WRITEC 6,*) ' ' 

PRINT ALL ARRAYS 
WRITEC6,'C2X, "ARRAY TSCA")') 

CALL ARRAYCMXMOS ,MXLOS ,MXGRD ,NMOS ,NLOS ,NGRD ,TSCA) 

WRITEC 6, 'C2X, "ARRAY TSCAM")') 

CALL ARRAYCMXMOS, MXLOS,MXGRD,NMOS,NLOS,NGRD, TSCAM) 

WRITEC6,'C2X, "ARRAY SHREB")') 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , SHREB ) 

WRITEC6,'C2X, "ARRAY SHRTEB' ')') 

CALL ARRAYCMXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , SHRTEB) 

WRITEC 6 ,' C 2X ," ARRAY PMLE ")' ) 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , PMLE ) 

WRITEC6, 'C2X," ARRAY PJ")') 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS ,NGRD ,PJ) 

WRITEC6,'C2X, "ARRAY PJM")') 



MAR08890 
MAR08900 
MAR08910 
MAR08920 
* MAR08930 
MAR08940 
MAR08950 
MAR08960 
MAR08970 
MAR08980 
MAR08990 
MAR09000 
MAR09010 
MAR09020 
MAR09030 
MAR09040 
MAR09050 
MAR09060 
MAR09070 
MAR09080 
MAR09090 
MAR09100 
MAR09110 
MAR09120 
MAR09 130 
MAR09140 
MAR09150 
MAR09160 
MAR09170 
MAR09180 
MAR09190 
MAR09200 
MAR09210 
MAR09220 
MAR09230 
MAR09240 
MAR09250 
MAR09260 
MAR09270 
MAR09280 
MAR09290 
MAR09300 
MAR09310 
MAR09320 
MAR09330 
MAR09340 
MAR09350 
MAR09360 
MAR09370 
MAR09380 
MAR09390 
MAR09400 
MAR09410 
MAR09420 
MAR09430 
MAR09440 
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CALL ARRAY( MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , P JM) 
WRITE(6,'(2X, "ARRAY PJl")') 

CALL ARRAY( MXMOS , MXLOS , MXGRD , NI-IOS , NLOS , NGRD , PJl ) 
WRITE( 6 , ' ( 2X , " ARRAY FJIM " ) ' ) 

CALL ARRAY( MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , P JIM) 
WRITE(6,'(2X, "ARRAY PEB")') 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , PEB ) 
WRITE ( 6, '(2X, "ARRAY PTEB")’) 

CALL ARRAY( MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , PTEB) 
WRITE(6, '(2X, "ARRAY PEBM")') 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , PEBM ) 
WRITE(6,'(2X,"ARRAY PTEBM")') 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , PTEBM) 
WRITEC6/C2X, "ARRAY V" ) ' ) 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , V) 
WRITEC6/C2X, "ARRAY XB")') 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , XB ) 
WRITEC6,’C2X, "ARRAY VT")’) 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , VT) 
WRITE C 6, 'C2X, "ARRAY XTB")’) 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , XTB ) 
WRITEC 6 , ' C 2X, " ARRAY XJ" ) ' ) 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , X J) 
WRITEC6/C2X, "ARRAY XJl")’) 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , XJl ) 
WRITEC6,'C2X, "ARRAY XEB")’) 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , XEB ) 
WRITEC6,'C2X, "ARRAY XTEB")') 

CALL ARRAYC MXMOS , MXLOS , MXGRD , NMOS , NLOS , NGRD , XTEB ) 

'k 

END 



SUBROUTINE WRMADC IYR,NYR,MU,KU,MO,KO,MD,KD) 

INTEGER IYRCNYR),KUCNYR),KOCNYR),KDCNYR) 

REAL MUCNYR) .MOCNYR) ,MDCNYR) 

WRITEC6,'C/9X,' 'DEVIATION' ' ,T26 , I2.T42, I2,T59 , 12/) ’ ) 

* CIYRCIW) ,IW=1,NYR) 

WRITEC 6, 'C9X, "UNDERAGE' ' ,T22,F9. 6, ' ' C ” ,12," )" ,T38,F9. 6, 

* " C",I2,")",T55,F9.6," C' ’ ,12,")")') 

* MUCD ,KUCD ,MUC2),KUC2) ,MUC3),KUC3) 

WRITEC 6, 'C9X, ' 'OVERAGE' ' ,T22,F9. 6," C " ,12,")" ,T38,F9. 6, 

* " C'M2,")",T55,F9.6," C' ’ ,12,")")') 

* MOCD ,KOCl) ,MOC2),KOC2) ,MOC3),KOC3) 

WRITEC6,'C9X, "TOTAL' ' ,T22,F9. 6," C" ,12, ")" ,T38,F9. 6, 

* " C" ,12,")" ,T55,F9. 6," CM2,")")') 

* MDC1),KDC1),MDC2),KDC2),MDC3),KDC3) 

WRITEC6,'C/)') 

RETURN 

END 



SUBROUTINE ARRAYC MM ,ML , MG , NM , NL , NG , OB J1 ) 

‘irkkkiHckkirkickkickkkickicirkkkkkkkickkificickkkkkic 

REAL 0BJ1CMM,ML,MG) 

DO 10 IG1=1,NG 



MAR09450 
MAR09460 
MAR09470 
MAR09480 
‘ MAR09490 
MAR09500 
MAR09510 
MAR09520 
MAR09530 
MAR09540 
MAR09550 
MAR09560 
MAR09570 
MAR09580 
MAR09590 
MAR09600 
MAR09610 
MAR09620 
MAR09630 
MAR09640 
MAR09650 
MAR09660 
MAR09670 
MAR09680 
MAR09690 
MAR09700 
MAR09710 
MAR09720 
MAR09730 
MAR09740 
MAR09750 
MAR09760 
MAR09770 
MAR09780 
MAR09790 
MAR09800 
MAR09810 
MAR09820 
MAR09830 
MAR09840 
^^AR09850 
MAR09860 
MAR09870 
MAR09880 
MAR09890 
MAR09900 
MAR09910 
MAR09920 
MAR09930 
MAR09940 
MAR09950 
MAR09960 
MAR09970 
MAR09980 
MAR09990 
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DO 10 IL1=1,NL 

10 WRITE(6/(2X,15(F7. 3,1X))' ) (0BJ1( IMl , ILl , IGl) , IM1=1 ,NM) 

WRITE(6, '(/)') 

RETURN 

END 



SUBROUTINE nT(CINV,Y, X) 

REAL CINV,X,Y 
TEMP=Y/( 1+CINV) 

TEMP1=(1.+Y)/(1+CINV) 

TEMP=SQRT(. 5+CINV)*. 5*( ASIN( -1. +2. *TEMP) + 

* ASIN(-'1.+2.*TEMP1)) 

X=TEMP 

END 

SUBROUTINE UOM( VINV , ACT , OBJ , OVER , UNDER , KO , KU , KD) 

^nHr*Tt************')Hr***^Wr****^V**iV^V*****^V****'AnVTV**^\r********VWr****Vr****** 



REAL OBJ 

IF (ACT .GT. OBJ) THEN 

UNDER = UNDER+VINV*(ACT-OBJ) 

KU=KU+1 

ELSE 

OVER = OVER +VINV*(OBJ-ACT) 

■ KO=KO+l 
END IF 
KD=KD+1 
RETURN 
END 

REAL FUNCTION RH0(XB ,XBB ,SSE, SSB ,DEE,KK,T,NN) 



INTEGER T 

REAL*8 SIGSQ, V, U 
SIGSQ=SSE/(NN+2) 

V=SSB/T 

U=KK*(KK-3)*(XB-XBB)**2/(V*(KK-1)) 
RHO=DMIN1(1.DO, DEE/DSQRT(U) ) 

END 



REAL FUNCTION MILLER( X, AVGINV ) 

TM = X/SQRT(AVGINV+.5)+l. 570796 
TEMPI = ASIN(SQRT(1. /(AVGINV+1))) 
IF (TM .LT. TEMPI) THEN 
MILLER = 0. 0 
GOTO 1 
ENDIF 

IF (TM .GT. 3. 14159-TEMPI) THEN 
MILLER = 1. 0 
GOTO 1 
ENDIF 

IF (COS(TM) .EQ. 0.0) THEN 
MILLER = 0. 5 
GOTO 1 
ENDIF 



MARIOOOO 
MARlOOlO 
MAR10020 
MAR10030 
MAR10040 
MAR10050 
MAR10060 
MAR10070 
MAR10080 
MAR10090 
MARlOlOO 
MARlOllO 
MAR10120 
MAR10130 
MAR 10 140 
MAR10150 
MAR10160 
MAR10170 
MAR10180 
MAR10190 
MAR10200 
MAR10210 
MAR10220 
MAR10230 
MAR10240 
MAR10250 
MAR10260 
MAR10270 
MAR10280 
MAR10290 
MAR10300 
MAR10310 
MAR10320 
MAR10330 
MAR10340 
MAR10350 
MAR10360 
MAR10370 
MAR10380 
MAR10390 
MAR10400 
MAR10410 
MAR10420 
MAR10430 
MAR 10440 
MAR10450 
MAR10460 
MAR10470 
MAR10480 
MAR10490 
MAR10500 
MAR10510 
MAR10520 
MAR 105 30 
MAR10540 
MAR10550 
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IF (COS(TM) .GT. 0.0) SIGNUM = 1. 

IF (COS(TM) .LT. 0.0) SIGNUM = -1. 

TEMP2=1-(SIN(TM)+(SIN(TM)-(1. /SIN(TM) ) ) /AVGINV)**2 
IF (TEMP2. LT. 0. 0) TEMP 2=0. 0 
MILLER=. 5*( 1. -SIGNUM* *TEMP2*’'f. 5) 

1 END 

REAL FUNCTION SCINV(X,A) 

-1 I-- ?-_>_■ t--t t- t. ■ t - - 1 - t I t -« ,t. — .t ..’ i-t- > ■ ■ -- t t- t-t . t - ■ It 1 < . f . t - 1. ■ 1. 1 ..t ti.t t t 

** SCALE INVERSION (TRANSFORM => ORIGINAL SCALE) 

R=X/SQRT(A+. 5) 

IF(R .LT. -1.570796) THEN 
SCINV=0. 

ELSE IF(R .GT. 1.570796) THEN 
SCINV=1. 

ELSE 

SCINV=. 5*(1. +SIN(R)) 

ENDIF 

END 



REAL FUNCTION VAR( ZZ ) 









DATA AA/1.6835/, B1/-.8934/, B2/.8991/ 
IF ( ZZ .GT. 2.2) THEN 
VAR = 1. 

ELSE 

IF ( ZZ . LT. 1. 001) THEN 
VAR = . 05 
ELSE 

VAR = AA*(ZZ**B1)*(ZZ-1. )**B2 
ENDIF 
ENDIF 
END 



REAL FUNCTION CHISQR( ID,0BJ,INV,ACT) 



*****Vf****V.-*Vf*******yfVf****Vr**********Vc********TWr****Vf*** ************** 



REAL OBJ,INV 

IF (ID. NE. 0 .AND. OBJ. NE. 0. .AND. OBJ. NE. 1. ) 

* CHISQR = (INV*(ACT-0BJ)**2)/(0BJ*(1. -OBJ)) 

END 

*******’>Wr*>V**********Vf***********yc*******************Vc *** * ************ 

SUBROUTINE READER(MOS, LOS, GRADE, YR,INV,Y, NMOS,SMOS, 

* NGRD,SGRD, SL0S1,SL0S2, ST1,ST2, IM,IL,IG,IT, lEOF) 

** READ RECORD AND ACCUMULATE LOSSES 

INTEGER MOS, LOS, GRADE, YR, INV, LOSS(8) 

REAL Y 

INTEGER STl, ST2 
INTEGER SLOSl, SLOS2 
INTEGER SMOS(40), NMOS 
INTEGER SGRD(IO), NGRD 
ICNT=0 

1 READ(1,100,END=999) YR,MOS ,GRADE, LOS , INV, (L0SS( I ) , 1=1 , 8) 

100 FORMAT(4I2,9I5) 

** CHECK IF RECORD MEETS SELECTION CRITERIA. OTHERWISE REJECT. 

IM=0 



MAR10560 
MARIOS 70 
MAR10580 
MAR10590 
MAR10600 
MAR10610 
MAR10620 
MAP.10630 
MAR10640 
MAR10650 
MAR10660 
MAR10670 
MAR10680 
MAR10690 
MAR10700 
MAR10710 
MAR10720 
MAR10730 
MAR10740 
MAR10750 
MAR10760 
MAR10770 
MAR 10 7 80 
MAR10790 
MAR10800 
MAR10810 
MAR10820 
MAR10830 
MAR 10 840 
MAR10850 
MAR10860 
MAR10870 
MAR10880 
MAR10890 
MAR10900 
MAR10910 
MAR10920 
MAR 109 30 
MAR10940 
MAR10950 
MAR10960 
MAR10970 
MAR10980 
MAR10990 
MARI 1000 
MARI 10 10 
MAR11020 
MAR11030 
MARI 1040 
MAR11050 
MARI 1060 
MARI 1070 
MAR11080 
MAR11090 
MARlllOO 
MARllllO 
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DO 10 I=1,NM0S 


MAR11120 




IF(M0S .EQ. SMOS(I)) THEN 


MAR11130 




IM=I 


MARI 1140 




GO TO 20 


MAR11150 




END IF 


'MARI 1160 


10 


CONTINUE 


MAR11170 




GO TO 1 


MAR 11 180 


*■ 




MAR11190 


20 


CONTINUE 


MAR11200 




DO 30 I=1,NGRD 


MAR11210 




IF( GRADE .EQ. SGRD(I)) THEN 


MAR11220 




IG=I 


MAR11230 




GO TO 40 


MAR11240 




END IF 


MAR11250 


30 


CONTINUE 


MAR11260 




GO TO 1 


MAR11270 


40 


CONTINUE 


MAR11280 


* 




MAR11290 




IFCLOS .LT. SLOSl .OR. LOS . GT. SLOS2) GO TO 1 


MAR11300 




IL-LOS-SLOSl+1 


MAR11310 


** 


YEARS OVER ST2 ARE USED FOR VALIDATION 


MAR11320 




IF(YR .LT. STl) GO TO 1 


MARI 1330 




IT=YR-ST1+1 


MARI 1340 


irk 


COMPUTE TOTAL LOSS 


MAR11350 




Y=0 


MAR11360 




DO 50 1=1,8 


MAR11370 




Y=Y + LOSS(I) 


MAR11380 


50 


CONTINUE 


MAR11390 


k 




MAR 11400 




RETURN 


MAR11410 


999 


IEOF=l 


MARI 1420 




END 


MAR11430 



\ 
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APPENDIX C 



SIMULATION FOR CHOICE OF AVERAGE 



INVENTORY VALUES 



PROGRAM MCSIM 






* SIMULATION COMPARISON OF INVENTORY MEAN PERFORMANCE * 

* * 



* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 

* 



PURPOSE: TO IDENTIFY THE BEST MEAN VALUE (ARITHMETIC, * 

GEOMETRIC AND HARMONIC) FOR USE IN THE INVERSION OF * 
TRANSFORMED DATA TO THE ORIGINAL SPACE- * * 

PROGRAMMER: CAPT C R DICKINSON, USMC * 

VARIABLE DESCRIPTION: * 

LAMBDA - PARAMETER USED IN THE POISSON RANDOM NUMBER * 
GENERATOR * 

TRIALS -* VECTOR OF INVENTORY VALUES GENERATED FROM * 
IHE POISSON RANDOM NUMBER GENERATOR * 

REPS - NUMBER OF INVENTORY VALUES TO BE GENERATED * 
BASIC - BASIC INVERSION OF TRANSFORMED VALUES * 

FTEl - FREEMAN-TUKEY EXACT INVERSION WITH * 

* 



^ FTE2 - FREEMAN-TUKEY EXACT INVERSION WITH 

* 



* 

* 



* 

* 

'k 

k 

k 

k 

k 

k 

k 

k 

k 

k 



Z - TRANSFORMED SCALE VARIABLE * 

MEANS - VECTOR OF ARITHMETIC, GEOMETRIC AND HARMONIC * 



MEANS * 

TRSUM - SUMMATION OF TRIALS FOR USE IN THE * 

ARITHMETIC MEAN * 

TRMULT - MULTIPLICATION OF TRIALS FOR USE IN THE * 

GEOMETRIC MEAN * 

TRINV - SUMMATION OF THE INVERSE OF TRIALS FOR USE * 
IN THE HARMONIC MEAN * 

INPUT/OUTPUT: * 

INPUT - HARD CODED PARAMETERS * 

OUTPUT - SIM OUT A1 (DISK FILE 03) * 



kkkkkkkkkk’i^rkkkickkkit'^ 



kkkkkkkki r kkkkkkkick kkkk kk kkkkkk kkk 



*** VARIABLE DECLARATION 

INTEGER TRIALS(IOO) , Y( 100) ,REPS( 3) ,PROB ,R, 

C COUNT, TRSUM, YTEMP.YSUM 

REAL Z( 100) ,T( 100) ,MEAN(3) ,BASIC(3) ,DIFF(3) ,LAMBDA(5) , 
C P(4),FTE1(3),FTE2(3), TRMULT 

DOUBLE PRECISION DPSEED 

* 

CALL EXCMS( 'FILEDEF 02 DISK REPLIC LISTING Al') 

* 

*** VARIABLE INITIALIZATION 

PI=3. 14159 
DPSEED=889246.D0 
LAM=5 

R=3 '■ 

PR0B=4 

COUNT=2 
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MCSOOOlO 

MCS00020 

MCS00030 

MCS00040 

MCS00050 

MCS00060 

MCS00070 

MCS00080 

MCS00090 

MCSOOlOO 

MCSOOllO 

MCS00120 

MCS00130 

MCS00140 

MCS00150 

MCS00160 

MCS00170 

MCS00180 

MCS00190 

MCS00200 

MCS00210 

MCS00220 

MCS00230 

MCS00240 

MCS00250 

MCS00260 

MCS00270 

MCS00280 

MCS00290 

MCS00300 

MCS00310 

MCS00320 

MCS00330 

MCS00340 

MCS00350 

MCS00360 

MCS00370 

MCS00380 

MCS00390 

MCS00400 

MCS00410 

MCS00420 

MCS00430 

MCS00440 

MCS00450 

MCS00460 

MCS00470 

MCS00480 

MCS00490 



DATA REPS/20,40,60/ 

DATA LAMBDA/1. ,2. ,4. ,8. ,16./ 

DATA P/.05,. 1,.2,.4/ 

* 

*** HEADER PRINT 

WRITE (02, ’(/T6, "LAMBDA' ' ,T15, "REPS' ' ,T27, "P' ' , 

C T34, ' 'METHOD' ' ,T46, ' 'ARITHMETIC' ' ,T59, 

G " GEOMETRIC", T71,‘' HARMONIC"/)’) 

* 

*** LOOP FOR LAMBDA PARAMETER 

DO 500 1=1, LAM 
DO 450 J=1,R 

CALL GGPOS ( LAMBDAC I ) , DPSEED , REPS ( J) , TRIALS , lER) 



10 

* 



c 

c 



* 



★ 

7CTC7T 



c 

c 



20 

* 



TrTTTr 



NSUM=0 

DO 10 J1=1,REPS(J) 

TRIALS(J1)=TRIALS(J1)+1 
NSUM=NSUM+TRI ALS( J1 ) 

LOOP FOR PROBABILITY PARAMETER 
DO 400 K=1,PR0B 
C0UNT=C0UNT+4 

LINE COUNTER FOR PAGE BREAK 
IF (COUNT. GT. 55) THEN 
WRITE (2,'(1H1)') 

COUNT=2 

WRITE (02, '(/T6, "LAMBDA’ ' ,T15, "REPS' ’ ,T27,"P' ' , 
T34,' 'METHOD' ' ,T46 'ARITHMETIC ’ ' ,T59, 
"GEOMETRIC' ’ ,T71, ' 'HARMONIC' '/) ' ) 

ENDIF 

TRSUM=0 

TRINV=0 

TRMULT=0. 

ARCSUM=0. 

TSUM=0. 

YSUM=0 

LOOP FOR REPITITION PARAMETER 
DO 20 K1=1,REPS(J) 

CALL GGBN(DPSEED,1,TRIALS(K1),P(K),YTEMP) 

YSUM=YSUM+YTEMP 

TRSUM=TRSUM+TR I ALS ( K 1 ) 

TRMULT=TRMULT+ALOG( REAL( TRIALS ( K1 ) ) ) 
TRINV=TRINV+1. /REAL(TRIALS(K1) ) 

Y(K1)=YTEMP 

Z(K1)=. 5*(TRIALS(K1)+. 5)**. 5*(ASIN(2. *Y(K1)/ 
(TRIALS(K1)+1. )-l. )+ASIN(2.*(Y(Kl) + l. )/ 
(TRIALS(K1)+1. )-l. )+PI) 
T(K1)=Z(K1)/((TRIALS(K1)+. 5)**. 5) 

ARC SUM=ARC SUM+Z ( K 1 ) 

TSUM=TSUM+T(K1) 

ZAVG=ARCSUM/REAL( REPS( J) ) 

TAVG=TSUM/REAL(REPS( J) ) 

MEAN VALUE CALCULATION 
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MCS00500 
MCS00510 
MCS00520 
MCS00530 
‘ MCS00540 
MCS00550 
MCS00560 
MCS00570 
MCS00580 
MCS00590 
MCS00600 
MCS00610 
MCS00620 
MCS00630 
MCS00640 
MCS00650 
MCS00660 
MCS00670 
MCS00680 
MCS00690 
MCS00700 
MCS00710 
MCS00720 
MCS00730 
MCS00740 
MCS00750 
MCS00760 
MCS00770 
MCS00780 
MCS00790 
MCS00800 
MCS00810 
MCS00820 
MCS00830 
MCS00840 
MCS00850 
MCS00860 
MCS00870 
MCS00880 
MCS00890 
MCS00900 
MCS00910 
MCS00920 
MCS00930 
MCS00940 
MCS00950 
MCS00960 
MCS00970 
MCS00980 
MCS00990 
MCSOIOOO 
MCSOlOlO 
MCS01020 
MCS01030 
MCS01040 
MCS01050 



* 









* 

30 



40 

* 

400 

450 

500 



MEAN( 1)=TRSUM/REAL(REPS( J) ) 

MEAN(2)=EXP(TRMULT’H1. /REAL(REPS( J) ) )) 

MEANC 3 )=REAL(REPS( J) ) /TRINV 

IF (COS(TAVG).GE. 0. 0) THEN 
SIGNUM=1. 

ELSE 

SIGNUM=-1. 

END IF 

DO 30 K2=l,3 

TAVG1=ZAVG/(MEAN(K2)+. 5)**. 5 
IF (COS(TAVGl).GE. 0. 0) THEN 
SIG1=1. 

ELSE 

SIG1=-1. 

ENDIF 

BASIC AND FTE INVERSE METHOD CALCULATION 

FTE1(K2)=. 5*( 1. -(SIGNUM*( 1. -( ABS(SIN(TAVG) + 

C (SIN(TAVG)-(l. /SIN(TAVG)))/MEAN(K2)))**2. 5)) 

FTE2(K2)=. 5*(1. -(SIG1*(1. -(ABS(SIN(TAVG1) + 

C (SIN(TAVGl)-(l. /SIN(TAVG1)))/MEAN(K2)))**2. )**. 5)) 

BASIC(K2)=. 5*(1.+SIN(2AVG/((MEAN(K2)+. 5)** 5) 

C -CPI/2. ))) 

OUTPUT STATEMENTS 
DO 40 K3=l,3 

IF (K3.EQ. 1) THEN 

WRITE (02,'(T6,F4. 1 ,T15 , 13,T25 ,F5. 3,T35 , " BASIC" , 

C T46,F7.4,T59,F7.4,T71,F7.4)') 

C LAMBDA(I) ,REPS( J) ,P(K) , (BASIC( IK) , IK=1 , 3) 

ELSE 

IF (K3. EQ. 2) THEN 

WRITE (02,'(T35,"FTE1" ,T46 ,F7. 4,T59 ,F7. 4, 

C T71,F7.4)') (FTE1(IK),IK=1,3) 

ELSE 

WRITE (02, '(T35,"FTE2’ ' ,T46,F7. 4,T59,F7.4, 

C T71.F7.4/)') (FTE2(IK) ,IK=1,3) 

ENDIF 

ENDIF 

CONTINUE 

CONTINUE 

CONTINUE 

CONTINUE 

STOP 

END 
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MCS01260 
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MCS01290 

MCS01300 
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MCS01330 
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MCS01420 
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